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@ A technique for the detection and removal of local defects in digital continuous-tone images. 

@ The present invention is a method for automatically detecting and correcting a wide range of local digital 
image defects with minimal user intervention. The detection process employs brightness and color thresholds in 
conjunction with magnitude thresholds on residuals of nonlinear spatial filters to separate defects from scene 
content with minimal confusion. The detected defects are then cosmetically corrected by combinations of 
nonlinear smoothing and grey-scale erosion. Several options are outlined for the feature selection, detection, and 
cleaning operations depending on source type and computational constraints. 
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Cross-reference to Related Applications 

The present application is related to U. S. Patent Application Serial No. 07/934.089 entitled "Process for 
Detecting and Mapping Dirt on The Surface of a Photographic Element" filed 08/21/92 by Robert Gray el 
5 al., corres. to EP 0 585 759 A1. 

Field of the Invention 

The present invention is related to the field of scanning images, primarily from photographic film, to 
10 form electrical equivalents of the scanned images and more particularly to the detection and the removal of 
defects due. for example, to dirt contamination or physical damage of the film. 

A portion of the disclosure of this patent document contains material which is subject to copyright 
protection. The copyright owner has no objection to the facsimile reproduction by any one of the patent 
disclosure, as it appears in the Patent and Trademark Office patent files or records, but othenwise reserves 
15 all copyright rights whatsoever. 

Background of the Invention 

Digital images created from electronic scanning of continuous-tone photographic film often reveal 

20 visually objectionable defects due to dirt contamination or surface damage of the film. These defects 
usually must be digitally corrected by manual digital retouching (e.g. via "dustbusting" and other cloning 
techniques), which require that the operator visually identify each local defect in the image. 

A patent of interest for its teaching in this art is U.S. Patent No, 4.189.235. entitled "Test Device For 
Dynamically Measuring The Degree of Dirt Accumulation On Bank-Notes" by Guter et al. That patent 

25 describes a method for inspecting opaque web materials for dirt accumulation. A reflected signal from the 
web material is sensed by three adjacent photosensors as the material moves past the sensors. When the 
signal sensed by the center sensor is significantly different than from that sensed by the side sensors the 
sensed area is classified as containing dirt. 

In conventional optical printing systems, if dirt or a scratch appears on the photographic image, the 

30 photographic negative or slide may be cleaned and the print remade or the print itself may be retouched. In 
the case where the images are written to an optical disc, it is not convenient nor desirable to rewrite the 
image, because many images are written onto a disc before the Images are viewed, and once an image is 
written it cannot be removed. Therefore, it is desirable to monitor the dirt present in photofinishing 
environments and to assess the effectiveness of any selected film cleaning methods. 

35 In U.S. Patent No. 4,907,156, entitled "Method and System for Enhancement and Detection of Abnormal 
Anatomic Regions in a Digital Image" by K. Doi et al. there is disclosed a technique for taking the 
difference between an "enhanced" and a "degraded" version of the input image to identify and remove 
structured anatomic background. A series of single thresholds are then applied to the difference images 
using the shape/size behavior of image regions to perform a feature extraction that is related to the 

40 detection of lung nodules. The feature extraction is implemented by the use of thresholds against which 
difference image pixels are compared to identify which super threshold "blobs" are lung nodules and which 
are not. 

A significant problem in the detection of dust is to provide an adequate way of detecting "small" 
objects. There are several possible approaches to this problem. For example, one could color-classify all 

45 the color pixels in an image, then spatially connect the pixels of a common color-class into segments, and 
then examine the spatial dimensions of the segment. This approach, while possible, has the drawback that 
the color classification must be quite accurate, or else some marginally-colored pixels will be misclassified. 
The problem is that differences in the direct RGB (or luminance-chrominance) image values do not alone 
give much contrast information between, say. a white small dust artifact and a bright grey extended scene 

50 region on which the artifact is superimposed. Thus after color-classification one may be left with very large 
segments that are reasonably neutral and bright, but which now mask visible dust artifacts. 

A better approach prior to classification is to preselect only those points which have large local contrast 
from their surround. This may be directly accomplished by (1) spatially smoothing the color image and (2) 
subtracting the smoothed color image from the original color image. The resulting difference (or "residuals") 

55 image will emphasize all color pixels which have a large magnitude compared to some weighted local 
average magnitude. The key decision here is then the choice of spatial smoother. The most obvious first 
candidate is a low-pass linear filter (e.g. a simple spatially weighted averaging of all grey levels (GLs) in a 
local window); the properties may be easily analyzed mathematically, and (perhaps more important) such 
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filters may be rapidly implemented in software and hardware (the residual image is basically a high-pass- 
filtered version of the original). Unfortunately, linear residuals images have the property of retaining ail local 
edges in an image, including very long edges due to scene discontinuities. The result is a very "cluttered" 
image (i.e. the residuals image contains many bright pixels having no connection with dirt or scratches); a 
significant amount of additional time must then be spent in scene reasoning to clean up the residuals map 
before defect decisions can be made. In short, linear-filter residuals poorly discriminate between local 
defects and true scene details. 

The identification of anomalies such as dust or scratches leads to the desirability of removing the same 
without the intervention of an operator. The present invention is directed to a methodology for identifying 
and cosmetically correcting such anomalies with minimal operator intervention. 

Summary of The Invention 

One of the preferred methods of the present invention is: A method for the detection and removal of 
defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of pixel 
values expected, of defects: 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range of 
expected defects: 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of residuals 
values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps b) 
and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

The present method, is called an autodustbuster (ADB) algorithm, and it has been successfully tested 
on both consumer-type color photographic images and color motion picture frames with good results. 

From the above it can ,be seen that it is a primary object of the present invention to provide an 
automatic local defect detection and removal technique for digital images. 

The above and other objects of the present invention will become more apparent when taken in 
conjunction with the following description and drawings wherein like characters indicate like parts and which 
drawings form a part of the present invention. 

Brief Description of The Drawings 

Figure 1 illustrates a system on which the present invention may be implemented. 
Figure 2 illustrates the windowing of an image raster to determine the grey level of a central pixel GLc 
at a position X. 

Figure 3 illustrates in chart form a histogram of GL values for an associated window. 
Figure 4 illustrates the averaging of the samples "nearest" in grey level GL to the central pixel GLc for 
the histogram of Figure 3. 

Figure 5 illustrates the replacement of the GLc value with the grey level average GLgv at the position X. 
Figure 6 illustrates the process flow for a multi-band original image. 
Figure 7 illustrates the process flow for a resolution pyramid of an image. 
Figure 8 illustrates the process flow as applied to an R, G, B image. 
Figure 9 illustrates a flow for ADB processing within a resolution pyramid. 

Detailed Description of the Invention 

Referring to Figure 1. the present invention is implemented using a digitizing scanner 10, a control and 
logic device 24, such as a computer with a display device 26. and an optical disc writer 36. As shown, a 
strip of film 15, containing one or more frames 18 of developed film, is placed into a scan gate 12 where 
the strip of film 15 is scanned under direction of control and logic circuitry 30. As each frame 18 is scanned 
the resultant scanned image data, represented by block 22, is digitized and transmitted to a memory 28 for 
storage. The computer 24 process the stored data, in a manner to be described, to provide output image 
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data 32 which may be written to an optical disc 34 by the optical disc writer 36 to provide a report as to the 
characteristics of the anomalies. The scanning device 10 is capable of quantizing pixel values into multiple 
brightness levels in separate red. green, and blue channels. A minimum number of brightness levels would 
be approximately 64 for an adequate quality image wilh a typical number being 256. Although the image 
data 22 has been described as originating from the operation of scanning film, it is well understood that 
other well-known techniques may be used to provide the image data. 

The autodustbuster (ADB) algorithm and several processing options that allow customization for 
different source types and computational constraints will now be described. 

The required ADB operations may be divided into (1) feature selection, (2) defect detection, and (3) 
defect correction. The operations must ensure that the maximum number of defects is detected and 
optimally corrected while minimally affecting scene elements which are not defects. Thus, the more 
confidence which the user has in the accuracy of the detection process, the more extreme can be the 
removal process; conversely, as one's confidence in avoiding "false alarms" from scene elements 
diminishes, milder removal methods may be used to avoid damaging real scene elements. 

Visible dirt artifacts in digital positives, created from scanned color negatives, possesses the following 
properties: they are "white" (low color saturation), "bright" (very high apparent exposure values), "small" 
(small spatial extent) and "contrast" (large differences from their immediate surround). (If dirty color 
positives were instead scanned, then the dirt artifacts would be "dark" rather than "bright"). In addition, 
many such fragments exhibit little brightness variation within individual defects. 

It is obvious that the above characteristics are rather loosely defined, and that many real scene 
characteristics may possess one or more of these attributes. Thus the goats of an ADB algorithm are to 
maximize the percentage of defects detected ("hits") while minimizing the number of real scene elements 
mistakenly considered to be defects ("false alarms"). These goals tend to be competitive, that is. setting 
thresholds on features to increase the percentage of "hits" tends to increase the number of "false alarms" 
and vice versa. A detection mechanism must be used in which defects appear as dissimilar as possible 
from true scene elements. 

In most imagery, spatially small objects which differ greatly in brightness or color from the immediate 
surround are relatively rare. The ADB algorithm exploits this property by the innovative use of a particular 
form of nonlinear spatial filter. The properties of this filter are described in the next section. 

Figure 2. illustrates in flow chart form the generic autodustbuster process flow. The original image of 
block 100 is a multi-band digital image represented by one or more matrices of pixel values or equivalentty 
grey levels (GL). In preprocessing block 110 some effective and convenient remapping of the original image 
may be performed. An example of preprocessing is the digital resampling of the image to reduce its size 
for computational expediency. An alternate preprocessing is the selection or recombination of the original 
image bands to produce a feature image with fewer bands. The detection suitability of the resulting feature 
image, in any case, should be at least comparable to that of the original image. Oval 120 represents the 
feature image that is outputted from block 110. The feature image is subjected to a band brightness test in 
block 130. The test applies low and high threshold evaluations to every pixel of each band in the feature 
image. Pixel values that are above the high threshold or below the low threshold in each band are 
considered to have passed the test and are forwarded to an edge-preserving-smoothing (EPS) filter 140. 
This filter non-linearly smoothes the selected pixel in all bands in a manner that retains major scene edges 
as subsequently described. The output of the filter at this pixel position is written to a corresponding 
position in the smoothed image 142. For those pixel values that are below the high threshold value, but 
above the low threshold value no subsequent EPS filtering is applied, i.e. the test is considered failed. Each 
band of the feature image pixel is copied to the resultant output of the smoothed image 142 and does not 
undergo EPS filtering. At difference node 145 the pixel values of each band of the smoothed image are 
subtracted from the pixel values at the corresponding position in the feature image 120 resulting in the 
residuals image 148. 

All of the pixel values from the residuals image are subjected to the residuals test of block 150. The 
residuals test consists of a minimum and a maximum threshold that is applied to each band of the residuals 
image 148. When bright defects are to be detected, for a pixel to pass all band values must exceed the 
minimum threshold and at least one band value must exceed the maximum threshold. When a dark defect 
is to be detected a pixel will pass if all band values are less than the maximum threshold and at least one 
band value is less than the minimum threshold. 

Pixels that pass are forwarded to a classmap 160 and marked with a distinctive value. The classmap 
160 is a single band image representation which is used to mark the location of detected defects. Pixels 
that fail are marked with a different value in the classmap 160. In a classmap conditioning block 170 the 
values within the classmap 160 may be further modified to accommodate known spatial characteristics of 
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the defects. For example, a one-pixel morphological dilation of the defect-marked (passed) pixels in the 
classmap is typically performed. The conditioned classmap 175 is directed to the cleaning block 180 long 
with the original image pixel values. Pixels in the original image that are not marked as passed in the 
classmap will be copied unchanged to the output cleaned image 190. Pixels which are marked as passed. 
5 I.e.. defect, are corrected using neighboring non-defect pixel values. The corrected values are then written 
to the cleaned image 190. 

A more detailed discussion of the major blocks of Figure 2 will now be undertaken. In ADB processing 
the feature image 120 used to detect the defect locations may be different from the original image 100 
which will be corrected. This is generally advantageous when computational time is at a premium and EPS 

10 filtering 140 of the full (generally three-color) image is too burdensome. In this case it is often possible to 
apply a preprocessing module 110 to create a one- or two-band "feature" image 120 which instead will be 
EPS-filtered and classified for defects. The effectiveness of this strategy will depend on the success with 
which a reduced-dimensionality feature image will retain the size, brightness, contrast, and color defect 
discrimination of the original image 100. This will of course be application-specific. For example, in digitized 

76 imagery originating from superimposed color cells (as in animated motion pictures) defects may be of any 
color and saturation; thus little is lost in detection performance using a reduced-dimension space which 
sacrifices color and saturation information. In this case good success has been achieved in using a "band- 
maximum" one-band feature image 120 band where every feature pixel value is the maximum value of the 
input red. green, and blue bands. By contrast, in Photo CD applications the principle source of defects is 

20 dirt on the negative during scanning, which appears as "white" specks in the digital positive; here saturation 
information should be retained as an important identifier, and two- or three-band EPS-filtered feature images 
are the norm. Figure 7 illustrates the retention of all three color original color bands for classification. 

Rescaling of the feature image grey-level range during preprocessing 110 is often advantageous. In 
current implementations this consists of linearly scaling to a 0-255 pixel value range even for 10- or 12-bit 

25 original data. This results in a slight speed improvement during EPS-filtering 140 when the histogram is 
searched at every window location (see Figures 4 and 5 along with associated discussions). In Photo CD 
applications (e.g. Figure 7 this rescaling may be determined by finding .5% and 99.5% cumulative 
histogram points of the color-corrected bands of the entire input image. The minimum .5%-point and 
maximum 99.5%-point from among all bands are mapped to 0 and 255 pixel values, respectively in the 

30 feature image, thus ensuring that no color shift occurs. If the input data is not color balanced it may be 
necessary to perform a color balancing step on the original data. 

EPS Filters and k-Nearest Neighbor (kNN) 

35 To detect the critical size of the desired defects a smoothing filter is desired which retains extended 
image edges but smoothes small local contrasty image regions (blobs). Such filters exist as various 
nonlinear noise-reduction methods under the category of "edge-preserving smoothers" (EPS). The most 
widely known EPS filter is the median filter, which is one member of the class of rank-order filters. In rank- 
order filters, the pixel magnitudes are sorted and an average of some contiguous number of the sorted 

40 magnitudes is computed. For example, in a median filter the "middle" element of the sorted magnitudes is 
retained (an average over one sample). The sample mean is, at the other extreme, also a rank-order filter 
(although a poor EPS); it provides the average of all sorted magnitude values. 

A poor property of standard rank-order filters for image operations is that they are insensitive to the 
relative spatial location of the pixels within the current window; e.g. the magnitude chosen as the "median" 

45 may exist at the center of the window or at an extreme edge of the window. A number of specialized 
nonlinear window operators have-been devised to account for local spatial structure. One example is the 
Sigma filter (J.-S. Lee. CVGIP.24,255-269. 1983). In this filter a histogram is created of the grey levels in a 
sliding window; based upon an estimate of the image noise, a range of histogram values about the grey 
level of the pixel at the center of the window is computed. One characteristic of the sigma filter, however, is 

50 that it is specifically designed for noise removal and its averaging range is solely dependent on the noise 
characteristics; the filter averages over a different number of pixels at each window location. For "blob- 
detection, however, what is desired is a fitter whose parameters are sensitive to local shape characteristics. 

A relative of the sigma filter which gives this control is the k-nearest-neighbor (kNN) filter (Davis and 
Rosenfeld, IEEE Trans. SMC, 7. 107-109 .1978). This filter is implemented in the EPS filter module 140 of 

55 Figure 2, This filter has the following operation per color band: 

Referring to Figure 3. slide a rectangular window of odd dimensions n x m over an image in raster 
order. At each location determine the grey level of the center pixel {GU) or equivalently pixel value, and 
also create a histogram of the grey levels in the local window as per Figure 4. Beginning at the histogram 
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bin containing the center grey level GLc, sum contents of the histogram bins on both sides of the starting 
bin in the following fashion: Test if the. contents in the bin for GU are greater than or equal to k pixels. Jf 
this is true, then write out as the averaged value GLgv just the existing value GLc. the contents for bin GL^ 
are less than k, then add the contents of the bin immediately below the bin and above GU to the sum and. 
5 again test whether the sum is greater than or equal to k. If the sum equals k, then calculate the average 
grey level according to the formula 

GLav = SUM {GLi X hist(i)} 

10 



SUM {hist (i) ) 



15 where the index i is over the summed histogram bins. If the sum is greater than k, then the extreme bins in 
the sum are proportionally reduced in weight in the calculation for GLgv so that the total bin weight in the 
sum approximately equals k (as shown in the attached source code). If the sum is less than k, then the next 
bin below GLc and above GU are summmed and tested against k. This process is repeated until the sum of 
histogram bins about GLc are greater than or equal to k. Continue until k pixels (k < n x m) have been 

20 summed. Find the average grey level of the summed pixels per Figure 5 and write this value to an output 
image at the location of the window center per Figure 6. 

The k-nearest neighbor filter has some interesting filtering properties which are not immediately obvious 
from the above description. First, any constant-GL image structure (blob) which has an area of k or more 
pixels in the local window will not be significantly smoothed: any structure of less than k pixels will be 

25 increasingly smoothed as k decreases. Second, at "ramp" edges between otherwise constant grey level 
regions the slope of the ramp will be increased (sharpened) by the filter. 

The first property is the one of significance for "blob" detection. For example, if the window size is 5 x 
5 pixels, and if k is, say, 15, then extended straight edges having a depth (constant pixel value in both 
directions perpendicular to the edge) greater than three pixels will be effectively unchanged, and thus 

30 disappear in a residuals image. "Blobs" of area less than 15 pixels that fit entirely in the window will be 
smoothed and thus will appear strongly in the residuals image. This is the desired residuals property for a 
local defect detector. 

A major concern for a k-nearest neighbor filtering is computational speed. This is not a linear filter and it 
is not separable, nor does it have other usable symmetry properties. The major operations per window 

35 location are histogram creation and then accumulation with testing of the total accumulated samples. The 
current k-nearest neighbor implementation does significantly reduce the cost of computing the histogram by 
merely updating the window histogram in a horizontal strip of image rather than performing complete 
calculations per window. This results in a fractional saving of approximately (n-2)/n in histogram-computing 
operations (for a n pixel x m line window), i.e. the savings increase as window length increases. (This fast 

40 histogram updating is equivalent to that used in the "fast-median" filter-Huang and Tang, IEEE Trans, 
ASSP. 27, 13-18, 1979). 

Residuals Test 

45 As mentioned above, a residuals image 148 is created as shown in Figure 2. The anomaly classifier 
steps contained in modules 150 and 170 accept as input the residuals image 148 and create as output the 
final classmap 175 which marks the pixels to be cleaned. 

The residuals test in module 150 creates an initial classmap by performing tests on each pixel of the 
residuals image 148. In the classmap conditioning module 170 a spatial segmentation (i.e. connected- 

50 component labeling) of the non-background pixels in the classmap may be optionally performed with 
additional tests applied to the same segment. 

In the simplest (single-state) version and for detecting bright defects, each color residual pixel is 
evaluated against two residual-GL thresholds. GLrmin and GL^max- m this implementation, all color residual 
pixels whose minimum band residual GL exceeds GL.^in andwhose maximum band residual GL exceeds 

55 GUmax a^o uniquely marked for subsequent cleaning; pixels which do not satisfy both of these constraints 
are marked to be left unaltered in the output image. (For single-band feature images GL^min equals GUmax)- 
The two-dimensional map of category markings per pixel is termed the "class map". These thresholds 
impose a minimum required contrast between a defect pixel and its immediate background. The (smaller) 

6 
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GLrmin threshold is generally set slightly above the "zero point" in the residuals image (e.g. to 133 for 
rescaled features with residuals biased by a value of 128 as described above) and acknowledges that valid 
bright defect pixels . e.g. due to dirt on a scanned negative, should not have a negative contrast (residual 
values less than the zero point) in any band. The (larger) GL^n^ax threshold states a requirement of 

5 significant local contrast in at least one band for the defect to be visible. (Recalling that these thresholds are 
on the residuals GLs. it is obvious that not all bands should be required to pass a large GUmax value; for 
example, a white dirt artifact against a bright red background would have large residual values in the green 
and blue bands, but only a small residual value in the red band.) 

For two-state anomaly classification, two sets of residuals thresholds are used, i.e. GUmini. GLnnaxi, 

/o GUmin2. and GU^axa- 'n general, the minimum-residuals thresholds are set equally, and slightly above the 
zero point as mentioned above; the maximum-residuals threshold is set somewhat higher for labeling high- 
confidence pixels than is the maximum-residuals threshold for the low-confidence pixels. This labeling 
results in a class map where the pixels possess one of three labels (normal, low-confidence dirt, and high- 
confidence dirt). In two-state classification, this class map undergoes image segmentation in module 170, 

75 i.e. all residuals pixels spatially contiguous in the class map which have the same class are. assigned a 
unique code value. (These code values can be thought of as uniquely labeling all the "blobs" of suspected 
defects in the image.) A fairly efficient algorithm for segmenting the class map is employed which involves 
finding the lowest prototype code value in a spatial "tree" of prototype segments of the same class, and 
then renaming all the segments in the tree with this lowest code value. 

20 Following segmentation, all segments which contain low-confidence defect pixels and which are 
spatially touching segments containing high-confidence pixels are upgraded to high-confidence (i.e. the 
class map values of the pixels are reassigned to the "high confidence of an defect" class). The purpose of 
this step is to make use of the adjacent presence of very obvious defects to decide that a slightly less 
contrasty region is part of the same defect. The segments which remain at low confidence are either 

25 relabeled as "clean" (i.e. all pixels of the segment are marked for no cleaning in the classmap) or else 
name their pixels marked for subsequent cleaning state in module 180. The segments, which are marked as 
high confidence of being defects, have their pixels marked in the classmap for subsequent erosion cleaning 
in module 180. 

30 Cleaning 

Two general cleaning methods are presently available for use within module 180 in the ADB algorithms: 
EPS substitution and greyscale erosion. In simple EPS substitution, the original input multi-band grey levels 
of image 100 are simply replaced by the band EPS-filtered grey levels of smoothed image 142 wherever a 

35 pixel is marked as a defect in the classmap 175 (as shown in Figure 2). This method has the advantage of 
computational simplicity in that EPS values have already been computed during defect detection; however It 
requires that the EPS image be stored for use during cleaning. In addition. EPS substitution is relatively 
forgiving of "false alarms", i.e. pixels classified as dirt artifacts but which are really part of the scene. 
Conversely, EPS-substitution cleaning may not completely remove all of an extended defect, to the degree 

40 that the defect itself influenced the local EPS values. 

In the greyscale erosion cleaning option, an iterative process occurs over detected defect regions, 
whereby GLs of defect pixels are replaced by a weighted average of original GLs of "normal" (i.e. non- 
defect or corrected defect) pixels within a small distance of the defect pixel. As the iterations progress, the 
effect of the "normal" GLs propagates toward the center of the original extended defect. The iterations end 

45 when all marked defect pixels have been replaced. Two precautions are in order. First, the "normal "-labeled 
pixels that are in immediate contact with the defect pixels at the beginning of the iterations are also 
relabeled to class "defect" prior to the onset of cleaning (via a morphological dilation operation which uses 
a 3x3-pixel probe function). This is because such pixels, while failing the detection thresholds, often have 
some contribution from the defect in their GLs which will be visible when propagated over multiple pixels 

50 during erosion. Cleaning their values as welt as those of "true" delect pixels can significantly improve the 
quality of the cleaning. The second precaution is that the spatial window over which the erosion weights are 
taken should be kept reasonably smalt, e.g. 5x5 samples. The reason for this is to minimize the influence 
of distant pixels which may not be part of the same scene structure as the region obscured by the central 
defect pixel. 

55 Current practice for the erosion cleaner is to use a 5x5-pixe! sliding window with a minimum of six 
clean-class contained pixels in order for a central defect-class pixel to be filtered; if the central pixel is 
defect-class, but fewer than six clean samples are contained, then the defect classification is maintained to 
the next iteration. If the central pixel is of defect class and six or more clean-class samples are contained in 

7 
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the window, then the GLs of the centra! pixel are replaced by a weighted average of the band GLs of the 
contained clean pixels, where the weights are linearly proportional to the inverse of the Euclidean spatial 
distance of the clean pixel fronn the central pixel. The cleaning is iterative (multiple passes through the class 
map) but is not recursive, i.e. the cleaning effect at one window location does not affect the cleaning at 

5 subsequent locations during the same iteration, in order to avoid directional biasing. 

Greyscale erosion has the advantage over EPS-substitution in being able to completely eliminate all 
trace of the artifact from the image; its disadvantages are increased complexity, a variable processing rate 
(the number of iterations is not known a priori), and the possibility for severe scene damage in the case of 
"false alarms" (i.e., real scene elements may be removed from the image). Despite these concerns, the 

10 excellent cleaning performance of erosion is currently the preferred method in both Photo CD (Figure 7) 
and Cineon (Figure 8) applications. 

Alternative Implementations: 

75 The ADB algorithm has several "flavors". These include (a) operation at a single spatial resolution, (b) 
operation within a resolution pyramid, (c) single-defect-state detection, (d) multiple-defect-state detection, 
and cleaning via (e) EPS substitution or (f) grey scale erosion. In the single-state versions, all pixels of the 
image are classified as either "non-defect" or "defect", i.e. only one defect state is treated. In multi-state 
versions, the "defect" classification is further divided Into subclassifications depending on the confidence 

20 that a pixel is indeed a defect; thus a three-state classification would consist of "non-defect", "high- 
confidence defect", and "transition-confidence defect". In single-state versions the cleaner used is either 
EPS-substitution or grey-scale erosion, but not both in one implementation. The multi-state version may 
employ a combination of both EPS-substitution for cleaning the lower-confidence defect regions and grey- 
scale erosion for cleaning the higher-confidence defect regions. 

25 The flavor combining options (a), (c). and (f) is a simple variant to implement, though it is not 
necessarily the fastest or best-performing version. In this case no spatial resolution pyramid is formed, the 
residuals image is not segmented, and the cleaner is chosen to be gray scale erosion for maximum 
correction. Figure 7 illustrates such a process flow and is typical of a Photo CD application. In this case the 
preprocessing 110 consists of rescaling the pixel values of the original image for maximum viewing contrast 

30 by linearly remapping the lowest .5-% band gray level to a value of zero and the highest 99.5% level to 
255. as previously described in the description for Figure 2. The brightness test 130 for dirt on scanned 
negatives consists of minimum threshold pixel values which must be equaled or exceeded in each band of 
the feature image 120 in order for the corresponding pixel to pass the test as a possible defect. The 
residuals classification 150 is a single state decision per pixel as either a non-defect or defect class, based 

35 upon a requirement that the smallest band residual equal or exceed the zero-background residuals level 
and that the largest band residual equal or exceed some user-determined threshold level which is above the 
zero background level. The classmap conditioning 170 in this case consists of a morphological dilation of 
the defect pixels in the classmap using a3x3or5x5 kernel. The cleaning 180 consists of greyscale 
erosion as previously described under the section called "cleaning." 

40 Figure 8 illustrates an alternate processing flow which retains the characteristics, (a), (c) and (f), is 
appropriate for removal of bright defects in Cineon applications. In this case it is desirable to minimize the 
processing time required by the large image format the processing is therefor directed at only detecting the 
largest defects and at minimizing the number of bands in the feature image. The preprocessing 110 in this 
instance consists of creating a reduced resolution image rescaling it to a zero to 255 greylevel range such 

45 that the color balance is corrected, and then extracting a single band feature image in which the output 
pixel value is the larger of the Red, Green and Blue input values. The residuals test 150 applies a single 
threshold which each residuals pixel must equal or exceed to be marked as defect. Classmap conditioning 
170 consists of morphological dilation as described above. The cleaning module 180 consists initially of 
performing greyscale erosion of the reduced resolution unrescaled input image as previously described. 

50 Next the cleaned reduced resolution multi-band image is bilineariy interpolated back to the original size and 
the conditioned classmap is replicated to the original image size. Finally those pixels in the original image 
which are marked "defect" in the replicated classmap are replaced by the pixel values in the interpolated 
cleaned reduced resolution image. 

It is possible to use ADB within a multi-resolution pyramid structure to provide some speed gains from 

55 use of multiple smaller EPS kernels. This implementation would also allow more effective cleaning from the 
use of EPS residuals thresholds which are tuned to the specific pyramid level. 

Figure 9 illustrates a flow for ADB processing within a resolution pyramid. In most respects, ADB 
cleaning within a spatial level is similar to that in a single-level ADB version, i.e. detection, classification, 
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and cleaning occur, and choices of one- or two-state and of EPS-substitution and/or greyscate erosion are 
the same. One key difference lies in the way tfie spatial pyrannid residuals innages are processed. 

In a spatial pyrannid. an image is represented generally as a low-resolution image plus a number of 
residuals images of increasing resolution. To reconslrucl the image to a particular resolution, one interpo- 

5 lates the low-resolution image by a prescribed method to the next resolution level and then adds the 
residuals values of that level to each pixel GL. This process continues until the desired resolution level has 
been reconstructed. The Photo CD file format is an example of a spatial pyramid structure. Note that these 
residual images are simply the pixel band difference between the original image at that resolution and the 
interpolated low-resolution version; they should not be confused with the EPS residual images discussed 

10 above. 

In a pyramid structure, all higher-resolution reconstructed pixels depend on the lower-resolution GLs at 
that spatial location. Thus when some subset of defect pixels is cleaned at a given resolution level, the 
residuals values at all higher resolutions at these locations must also be adjusted; failure to correct the 
residuals' values could result in "ghost" images of the original defects upon pyramid reconstruction. 

15 Although various residuals processing schemes are possible, the method currently used for ADB is the 
following relatively simple one: the base (lowest level) version of the original image is cleaned as in Figure 
7. The classmap 175 prior to base cleaning is then replicated in module 200 to the size of the next highest 
resolution level. The defect-class pixels in the replicated classmap are then dilated in module 210 by a 
spatial extent equal to the region of support of the pyramid smoothing function for that level, i.e., the half- 

20 width of the kernal used to create the next-lower-resolution level. The resulting dilated classmap 220 is then 
compared with the pyramid-residuals image 230; the pyramid residuals image is assigned a zero pixel 
value in every band wherever the corresponding pixel in the class map 220 is labeled as defect. The 
resulting modified pyramid-residuals image 250 is then added per-pixel in module 270 to a version of the 
base level cleaned image 190 which has been spatially interpolated to the present pyramid resolution level 

25 in module 260. 

This results in a reassembled (base + 1)-level image 100' in which the results of the base-level ADB 
learning have been incorporated. This image is then ADB-processed starting with the preprocessing step 
110', and the process of Figure 9 thus repeats until all levels have been reassembled and cleaned. Note 
that ADB cleaning may be omitted from higher-resolution levels at the expense of small defects not being 

30 corrected, but that all higher-level, pyramid-residuals must be zeroed as described in order to ensure that 
no defect "ghosts" will appear in the final image- 
While there has been shown what are considered to be the preferred embodiments of the invention, it 
will be manifest that many changes and modifications may be made therein without departing from the 
essential spirit of the invention. It is intended, therefore, in the annexed claims, to cover all such changes 

35 and modifications as may fall within the true scope of the invention. 
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Scanner device 
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Scan gate 
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Film strip 
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Film cleaner 
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Frame 
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Image data 


45 


24 


Control and logic device (computer) 




26 


Display device 




28 


Memory 




30 


Logic circuitry 




32 


Output image data 


50 


34 


Optical disc 




36 


Disc writer 




100 


Original image 




110 


Preprocessing 




120 


Feature image 
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130 


Band brightness test 




140 


Edge-preserving-smoothing (EPS) filter 




142 


Smoothed image 




145 


Difference or subtracting node 
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148 


Residual image 


150 


Residuals test 


160 


Ctassmap 


170 


Conditioning block 


175 


Conditioned classmap 


180 


Cleaning block 


190 


Cleaned image 


200 


Classmap replication block 


210 


Classmap dilation block 


220 


Replicated dilated classmap 


230 


Higher-resolution spatial resolution-residuals image 


240 


Spatial residuals zeroing block 


250 


Modified higher-resolution spatial resolution-residuals imago 


260 


Cleaned-image interpolation block 


270 


Pixel value summing node 


100* 


Rebuilt modified higher-resolution image 


110' 


Higher-resolution preprocessing 



The Invention may be summarized as follows: 

1. A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

2. A method for the detection and removal of defects in digital images comprising the steps of; 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

d) forming a residuals image which consists of the difference between the feature image pixel values 
and the filtered values of step c); 

e) testing each pixel value of the residuals image of step d) to determine if the residual value is within 
a range of residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

3. A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
residuals values expected of defects; 
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f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

4. The method wherein the map formed in step f) is conditioned by morphological dilation of the pixels 
marked as defect. 

5. The method wherein the map formed in step f) is conditioned by morphological erosion of the pixels 
not marked as defect. 

6. The method wherein the map of step f) is formed with each pixel marked according to the magnitude 
of its values according to steps b) and e). 

7. The method wherein the map formed in step f) is formed with each pixel marked according to the 
magnitude of its values according to steps b) and e)and is conditioned by morphological dilation of the 
pixels marked as defect. 

8- The method wherein the map formed in step f) is formed with each pixel marked according to the 
magnitude of its values according to steps b) and e) and is conditioned by morphological erosion of the 
pixels not marked as defect. 

9. The method wherein the filtering of step c) is performed with a k-nearest neighbor filter , 

10. A method for the detection and removal of defects in digital images represented by multiple spatial 
resolutions ordered from lower to higher, comprising the steps of: 

a) preprocessing a lower spatial resolution representation of the digital image by forming a feature 
image in which defect pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; 

g) correcting the lower resolution representation of the digital image using the map created in step f); 
and 

h) correcting the higher resolution representations of the digital image using the map created in step 
f)- 

11. The method wherein one or more of the higher resolution representations previously corrected using 
the lower resolution defect map of step f) are again processed using steps a) through h) to detect and 
correct for defects too small to be detected at lower spatial resolution representations. 

12. The method wherein the filtering of step c) is performed with a k-nearest neighbor filter. 

13. The method wherein the filtering of step c) is performed with a k-nearest neighbor filter. 

14. An apparatus for the detection and removal of defects in digital images comprising: 

means for preprocessing a digital image represented by pixel values to form a feature image in 
which defect pixels have large local contrast from their neighboring non-defect pixels; 

first means for testing the value of each feature image pixel value to determine if the value is within 
a range of pixel values expected of defects; 

filter means for edge-preserving spatial filtering of each of the feature image pixel values; 

means for forming a residual value for each pixel as a function of the difference between the 
corresponding feature image pixel value and the filtered value from said filter means; 

second means tor testing each residual value from said forming means to determine if the residual 
value is within a range of residuals values expected of detects; 

means for forming a map in which each pixel that is within the range of expected defects according 
to said first and said second testing means is marked as a defect; and 

means for correcting the digital image as a function of the formed map. 

15. An apparatus for the detection and removal of defects in digital images represented by multiple 
spatial resolutions ordered from lower to higher, comprising: 

means for preprocessing a lower spatial resolution representation of the digital image by forming a 
feature image in which defect pixels have large local contrast from their neighboring non-defect pixels; 

first testing means for testing the value of each feature image pixel value to determine if the value is 
within a range of pixel values expected of defects; 
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filter means for edge-preserving spatial filtering of each of the feature image pixel values that are 
within the range of expected defects; 

residual means for forming a residual value for each pixel as a function of the difference between the 
corresponding feature image pixel value and the filtered value of said filter means; 

second testing means for testing each residual value from said residual means to determine if the 
residual value is within a range of residuals values expected of defects; 

mapping means for forming a map in which each pixel that is within the range of expected defects 
from said first and said second testing means is marked as a defect; 

means for correcting the lower resolution representation of the digital image using the map created 
by said mapping means; and 

means for correcting the higher resolution representations of the digital image using the map created 
by said mapping means. 

16. A method for the detection and removal of local defects in digital images comprising the steps of: 

a) generating a feature image from an original digital image; 

b) creating an EPS-residuals image using the feature image; 

c) testing and mapping image pixels based upon the values of the EPS residuals image and the pixel 
values of the feature image; 

d) modifying the mapped image as a function of mapped values of neighboring pixels; 

e) cleaning of defect pixels by an EPS-substitution process or by a greyscale erosion process; and 

f) cleaning of defect pixels in different spatial resolution representations of the original image if 
different spatial resolution representations are present. 
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Busteps 042193: (4-21-93) Applies a selective-averaging 
(Tc-NN smoothing) filter to a 3-band image; also 

calculates the difference (residual image) between the original 
£ the smoothed images. 

Search window size must be symmetric, odd of dimensions 

3X3 thru 51x51. * ^ ^r.ese^ * 

Filter currently treats data as in range 0-4 999 *only . 
Input image array is stored as Int*2. 

Intended for operation on a Sun wor)cstation (e.g. Sparcstation 10) 
running Sun OS- 



Calls: 



20 



25 



30 



35 



40 



45 



50 



10 



Read Hdr isl . . , ^ ^ 

~ Timage-f ormat-specif ic calls) 

Ihister ^ 
(image-f ormat-specif ic ca^ls) 

Ihismnmx 
Winbust 042193 

~ (image-format-specific calls) 

Hist Init 
■Hist3Jpdate 
Iselavh 
Add Image 

" (image-f ormat-specif ic calls) 

program Busteps_042193 
integer*2 npix,nlin 

integer'*2 ibuff (4 0 96, 3) , buff n (4 096, 51, 3) 
integer nav, npp, diagf Ig 

Integer npixm,nlm, np,nl,band, nbands, nbx 
integer pixf ormi, pixf ormo, pixf ormro, dumml, dumm2 
integer gmin (3) , gmax (3) 
integer hist (5000, 3) 
integer*4 britelo, britehi 

real*4 f actl, f act2, bias 

real*4 deltime, Dtirae, timearray (2) 

character*50 inf il, outf il, resf il 

data npm. nlra/51, 51/ ! Max EPS window size! ^ w ^ i 

Tdata pix?or;ro/3/ ! residuals image -ritten as unsigned byte! 

data factl,fact2/1.0;-1.0/ !Used m Add_Image! 

data npixm,nbx/4096, 3/ 

data diagflg/0/ !For winbust3 ! 

common /par/dummJ., dumm2,nav 

Interactive parameter input 

write (6,5) 

formate "Busteps 042193- /, ^ -.k.^h*.'/ 

Tnput image may be i,2, or 3 bands;/. 
Input values should be in range 0-4999; /r 
Max. image width: 4096 pixels per line; /, 
Max. window height: 51 pixels'/. 
Left- & right -edge virtual mirroring /) 

write(6,10) ^ ^^^^ , 

format ( ' Enter in filename [ch50] 5) 

read(5> (a50)')infil !Must be ISL IPP imagefile! 

write (6,30) 



55 
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2 

^ T-P*; f ilenaine[ch501 : ' #5) 
read(5,' (aSO)') resrii /„dd >=3,-<51) . S) 

r=sti^"»"' -YisirsyKs 

Call Exitd) 



35 
40 



60 



70 



Endif 

nl=np 
npp=np*nl 



write (6, 50) averaging points in window 

write (d, > ' 



stop 
Endif 



-^-'V-'^J^ter low- . high-brig.tness threshold. [intl:',5, 
SIlNi^Sritelo.britehi 

^'"^1'- Inter out residuals bias constant: S) 
-^"'*riirisl(infil.npi.,nlin,nba.ds,pixfor..) 

" --!l3e i?'' ^nrut * too large-abort. 

Stop 



Endif 



initialize arrays s constants . 

npixin=npix 

„bx=nbands ^^^^^^ ^^^^^^l 

deltime=Dtinie (timearray) 

pixfonno=pixformi -Added 8/-;/ 92! 

pixfonnro=pixformi 

™av of input image 
Calc histogram ain.max oi i f 

write(6,l00)infil ^q, 
100 formate input file is . . ^ .. „ properly fil"r the image. 
Note: This min/max call is r qu 

150 format (/ 

DO ib-l^nbands 
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10 



15 



20 



35 



vrite(6,151) ib, grain (ib) , gmax (ib) 
151 f ormat (52t, i4,2i9) 

Enddo 

:.,tesc if image mins & maxs are within LOT range of 0-4999: 
Do ib=l,nbands 

If (gmax(ib) .gt. 4999) then 

write (6,*)' Warning— gmax ,gr. 4 999 STOP 

Call Exit(l) 
Elseif (gTBin(tb) .It. 0) then 

write(6,*)' Warning—gmin -It. 0 STOP' 

Call Exit(l) 

Endif 
Enddo 

Note- If gmax gt 4999, then either rescale input data, or redimension 
his arrays in selav3 & selavb, & reset threshold on above tester. 



Create EPS-filtered image 



write (6,*)' [Beginning filtering}' 
If (np .ge. 3 .and, np .le. npm) then 

Call Winbust 042193 (infil, outfil, ibuff, buffn, 

1 npixm, nEx, nlm, np, nl, gmin, graax, pixf ormo, britelo, briteni, 

2 diagflg) 

Else 

Write (6,*) ' Np out of range stop' 

Endif 



Create EPS residuals image 

25 Call Add Image (infil,outfiI,resfil,factl,fact2, bias, pixformro) 

Treslduals file is unsigned byte w/ mean="bxas ! 

write (6/ 200) outf il 
200 formate Ending; EPS file is ',a50) 

write (6, 201) resfil ^ , 

formate EPS residuals file is ',a50) 



201 

30 



deltime=Dtime (timearray) 
write (6, 800) timearray (1) , timearray (2) 
80C format (/' User time (sec):', 

X f8-2,'; system time (sec) : ' , f 8 . 2/) 

end 



Ibismnmx: (5-28-91) Calcs. min & max values in a multiband histogram. 

Subroutine Ihismnmx (hist , nbands, gmin, gmax) 
40 Integer hist (5000, 3) ,nbands,gmin (3) , gmax (3) 

Do ib=l,nbands 

Do Jc=l,5000 

If (hist()c,ib) .gt. 0) then 
gmin(ib)=)c -1 
go to 50 

45 Endif 
Enddo 

50 Do )c-5000, 1,-1 

If (hist()c,ib) .gt. 0) then 
gmax(ib)«)c -1 
go to 60 

50 



Endif 

Enddo 

60 continue 
55 Enddo 
return 
end 
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Read Hdr ISL: (9-25-90) Reads key info from header of an ISL-frat 
- image header. Closes image ^ ^l^_^F°^f f ^*^^l^f ^°2==r=-=== 

SuJDroutine Read_Hdr_ISL (filename, npix, nlin, nbd, pxf orm) 

Integer imgnum, ios , npx , nln, nbd, pxf orm, rwmode 
Integer * 4 opnimg , getdef 
Integer*2 npix,nlin 
.Character*50 filename 

I — — ' 

include ' /local/include/iopackagef ' 

rvmiode=0 
iingnum=0 

ios=Opnimg (imgnuin, filename, rwmode, - FALSE. ) 
If (ios -NE. SYSNRM) then ^ ^ 

write <€,*) ' Error during OPNIMG — Return 

return 

Endif 

ios=(^tdef (imgnuin, npx, nln, nt>d,pxrorm) 
If (ios .NE, SYSNRM) then 

write (€,*)' Error during Getdef — Return 

return 

Endif 

Call Clsimg (imgnura) 
npix=npx 
nlin=nln 

return 
end 
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- Ihister: (5-28-91) Calcs- multiband histogram of an image using a line 
: buffer only. Intended for use w/ seiectave2 ,f or . 

Subroutine Ihister (inf il, npixm, nbx, ibuf f ,hist) 
Character*50 infil 

integer npixm, nbx , . , . . 

integer arrtyp, ©status, gstatus, rstatus, iingin,npx,nln,nbands,pixfonn 
integer opnimg, getdef , rdline, line 
integer hist (5000, 3) 
integer*2 ibuff (npixm, nbx) 

include ' /local/include/iopac)cage . inc' 
arrtyp-IDINT2 

:.. .Initialize the histogram: 
Do ib-l,nbx 

Do )c=l, 5000 

hist{k,ib)=0 

Enddo 

Enddo 

: . - . en the image : 
imgin«0 

ostatus=opnimg (imgin, infil / 0 , .false . ) 
If (ostatus .ne. srSNKM) then 

■ write (6, 15) infil ^ ^ 

15 formate Ihister; input image ',a50,' not open — abort ) 

stop 

Endif 

gstatus=getdef (imgin, npx, nln, nbands^pixform) 
If (gstatus -ne- SYSNRM) then 

write(6,20) . ^ ^ 

20 formate Ihister: Gstatus error on input — abort ) 

stop 

Endif 

:,,.Read the image lines into the histogram: 
Do j=l,nln 

line^j -1 . _^ . 

rstatus=rdline (imgin, line, -1 , ibuf f , npixm, arrtyp) 
If (rstatus -ne. SYSNRM) then 

write <6, 50) j , . , . , , w » \ 

50 formate Ihister: error reading line ,i4, —abort ) 

stop 

Endif 

Do ib=l,nbands 

Do i=l,npx 

Jc=ibuff (i,ib) +1 

hist ()c,ib)=hist (k,ib) +1 

Enddo 

Enddo 

Enddo 

Call clsimg (imgin) 

return 

end 
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Minbust 042193: , H=.nHc: 

Uses a barrel-shifting line buffer in multiple bands. 
Implements a "fast-histogram" method. 



10 



: Key Parameters^^^^^ input image file; not overwritten. 

Ch*50, name of output filtered image tile. 
Int*2; all-band line buffer array; must be dimensioned by 

calling program, . . . . 

Int*2, strip buffer array; must be dimensioned by calling 

program. . , ^ . 

Inf4, Max. allowed width of image, m pixels per line. 

IntM! Actual width of image, in pixels per line. 

Read from image header. ^ , 

Inf4, Max. allowed no. of bands in input image. I.-*pUu. 
Int*4, Max. allowed height of filter window (in lines) 
& thus max, allowed height of strip buffer. 

IntM, Actual 'width, height of filter window in pixels, 
lines. Input. ^ . , i.. ^ r-r ^ 

gmin,gmax:lnf4, 3-element vectors. Assumed mm & °^^^^,^^^h inntir 
^ of input image. Assumed previously computed, input 

Note this limits current max. no. of bands to 3. 
pixformprlntM, Datatype of output filtered IPP ^^f |f ^^^^^A.^ 
E.g. 3=uns byte, 4-signed Int*2, 6-signea k 4. 
Input . 

amn,qmx: Int*4, Retur.ned via Common. Set to current 
^ -band value of gmin & gmax, respectively. These values 

be set prior to calling Iselavh; Function Iselavh uses 
these values passed via Common.. , , 

nav: $n"4, no of samples to average in each filter window. 
Passed thru Common by caller. 
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outfil: 
itouff : 

buffn: 

npixm: 

npix: 

nbx: 
nlm: 

np, nl : 



MUST 



General 
(1) 



(2) 



comments: . ^ ^. . - 

99% of the complexity of this subroutine is due to 

the implementation of the strip buffer. T^e„a^^^"^^?^.°5 

strip buffer is that it allows -convolution"-type window 

operations while keeping at any ^ one ^i^^ ^^T^fV^^^^y f__ 

npix*nl*nbands worth of image, instead of npix*nlin*nbands tor 

nS?^ tha^^'tEf E^riifter is *NOT* recursive, ile. it only uses 
K ortgina! Il values in its input weights; , thus output f^^^^ 
previous window filter do_not_aff ect_neighboring^wm 



Subroutine winbust 042193 (infil, outfil, ibuff, buff n, npixm, nbx, 
nim,np,nl,^in,gmax,pixformo,britelo,britehi,diagflg) 

Ctxaracter*50 infil, outfil 

integer imgin, irogout , ostatus , gstatus 

integer nbx,nbands ,pixf orm, pixf ormo „„i v noixb 

integer pstatus, ostatus2, rstatus, wstatus, npx, nln, npixm, npix, npixc 

integer nlin,nlinb 

integer arrtyp, lininc, jtype 

integer gmn, gmx, Iselavh, nav, gmm (3) , gmax (J) 

integer nlm, np, nl , nl2 , np2, nbinx, gl_cent 

integer opnimg, getdef , putdef , rdline, wrlme 

integer d-iagflg 

integer*4 britelo, br it ehi, count 

integer*2 ibuff (npixm, nbx) , buff n (npixm, nim, , nbx) 
lnteger*2 his (5000) 



50 



55 



22 



EP 0 624 848 A2 



11 



35 



lnt€ger*2 els (4096) !Note implicit max. of npix-4096 ! 

data nbinx/SOOO/ 
common /par/cron, gmx, nav 

np2-np/2 ! Half -width 



count=0 

rT'To^en'InputToutput images for sequential I/O: 



70 



80 



1 



imgin=0 

ostatus=sySIiRM , . , 
ostatus=opnxing(iingin, infil,0, .xaise.; 
If (ostatus .ne. SYSNKM) then 

?^^I^^'SUSst:042193:X^ image ',a50,' not open-abort') 

call Exit(l) 
Endlf 

|itttui:geS (imgin , npx, nln, nbands , pixf orm) 

If (gstatus .ne. SYSNRM) then 

write (6, 75) i„^„*- Anorr' ) 

formate Gstatus error on input aoorc ; 

Call Exitd) 

Endif 

npix-npx ! Actual no. pixels/line in image 

" '°'""5,:iS;S?5}i"!lSS=t_042193= »pl=.l. > «.!».-">■>"' 

Call Exitd) 

Endif 

nlinb=nlin +2*nl2 . , 

!nbands=Actual no. of bands m image 

?Smt^t/'°-Winbust 042193- : Npix Nlin Type Nbands Typeo Np NaV , 
;rf^l'j^?8S'iSix?n!!e^ii.for«nbands,pixfor«o,np,nav, 



81 



40 



1 britelo,brxtehi.diagflg 

format (21x, 2i5, i6, 17, i6, 15,217, iB) 

iingout=0 
pstatus=^SYSNRM ^ ^ 

pstatus=|utdef (imgout, npx, nln, nbands, pxxformo) 

(pstatus^.ne.^sySNHM,^th^^^^^ output-abort' 

Call Exit (1) 

45 Endif 

OStatus2-SYSNRM ^-F.i T fAlse ) 

ostatus2=opnimg SJ^^^^^^^he-inlil 

{ofitatus2 -ne. SYSNRM) then ,k«^' 

If <°«"^^f^g(6^,). Output image not open-abort 
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Endif 



12 

Call Exit(l) 

lininc=nlin/8 'Variable used in writing current line-number 
jtype-lininc ! Ditto (just screen i/o related) 
rstatus-SYSKRM 
wstatus=SYSNRM 
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: Load Ist strip into barrel buffer: 

Do i=l,nl2 •+! !Do per line in 1st buffer strip: 

line=j -1 ! (Damned -zero-starting* C array convention) 

rstatus=rdiine (imgin, line, -1/ ibuf f , npixia, arrtyp) 
If (rstatus .ne. SYSNRM) then 

write (6, 110) j , , c 

110 " formate "Winbust_0421 93" : Read error at line '-abort ) 

Call Exit(l) 
Endif 

! Transfer GLs from line buffer to strip buffer: 

! (Someone should be able speed up all this moving around of data) 
!i =horiz pixel index for original image; . , , ^ • w 

!i2=horiz pixel index for central (non-mirror) pixels in strip m 
!i3==horiz pixel index for edge (mirror) pixels in strip buffer. 

j2=j +nl2 
j3=nl2 +2 -j 
Do ib=l,nbands 

Do i==l,np2 +1 

i2=i +np2 

i3=np2 +2 -i 

buffn(i2, j2,ib) 

buffn(i3, j2,ib) 

buffn{i2, j3,ib) 

buffnti3, j3,ib) 

Enddo 

Do i-np2 +2rnpix -np2 - 
i2=i +np2 
buffn(i2, j2, ib)^ 
buffn(i2, j3,ib)^ 

Enddo 

. Do i=npix -np2 , npix 
i2=i ■*'np2 
i3=2*npix +np2 
buffn(i2, j2,ib) 
buffn(i3, j2,ib) 
buffn(i2, j3, ib) 
buffn (i3, j3, ib) 

Enddo 



ILeft 



=ibuff (i,ib) 
=ibuff (i,ib) 
=ibuff (i,ib) 
^ibuff (i,ib) 



J Central 



:ibuff (i,ib) 
:ibuff (i,ib) 

•Right 



=ibuf f (i, ib) 
=ibuff (i,ib) 
=ibuff (i,ib) 
=ibuff (i,ib) 



Enddo 

!Wrlte a status display to standard output! 

If ((j -nl2 -1) -ge. jtype) then w «i-Hin 

write (6, 600) j !Write current line number to staio . 

jtype=jtype +lininc 

Endif 

Enddo !End top "j" loop! 

jend=nl l^jend" is barrel index to most recent line 
jc«jend -ni2 !*jc" is barrel index to "middle" line 
Above sections only serves to fill the strip buffer for the 

first time; it gets more complicated as the buffer contents move 
"down" in the image. 
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13 



the buffer line index is rotated. 

jo=j -nl2 -1 linage andex of output line. 

iMark which pixels in the current line will be filtered: 
Do i^npa +l,npixb -np2 

if (buffnd, jcib) .GT. britelo .and filter' 
10 , buffn(i, jcib) .LT. britehi) then IDon t tiiter. 

^ cls(i)=0 

go to 210 

Endif 

Enddo , . , . , 

f-l<;H\=l !raark for filtering! , . , ^ 

louncUount *1 -Increment to-be-filtered pixel counter! 
15 210 continue 
Enddo 

Do ib=lfnbands 

gmn=gmin (ib) 
gmx^gmax (ib) 
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!Load window in central part ff f £!f 1 



i i^r^urrent band 
; Nit2°tht"^Ixe!I-within 1/2-wi.dow width of vert image boundaries 
iBeSow 'ii^ll^e-as above "12", i.e. strip buffer horiz index. 

- . !Center-of -window index, 1st window! 

io=i -np2 -output pixel index! 

Call Hist init (nlin,np,nl,np2,buffn-,npixn,nbx, 
npixb,nbands, i, jend, ib.his.nbmx) 

gl_cent:=buffn(i,jc.ib) = f ^f.-i^'^locatSon ! 



" "='=*i;ufnio'ibl=1selavh(his.nbinx,gl_cent) -Filter 
Else 

ibuff (io. ib) =gl_cent 

35 Endif 

[Update histogram for later window locations in a line: 

DO i=np2 +2,npixb -np2 !Center-of -window index! 
io=i -np2 'Output pixel index 1 

Call Hist update <nlin,np,nl,np2, buffn, " 
1 npixH, nbx,npixb,nbands, 1, Dend, ib, his^nbinx) 



!If passes brightness thresholds ^est, then filt« 
gl_cLt^buf fn (!, jc, ib) l^^^of^window^center 

^"i^^^^(inb}ilse!2vh(his,n^^^^ '.Filter 



Else 
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ibuff (io, ib)=gl_cent 

Endif 
Enddo 



Enddo !End band loop! 

line=jo -1 

wstatus=wrline (lingout, line, -1 , ibuf f , npixm, arrtyp) 
If (wstatus -ne. SYSNRM) then 

write (6, 410) jo ...... w , . 

formate -Winbust_042l93" : write error at line i4, ' -abort ' ) 

Call Exitd) 

Endif 

! Rotate barrel-buffer line indices: 

jend=jend +1 

If (jend .gt. nl) jend=l 

jc^jc +1 

If <jc .gt. nl) jc=l 

If <jo -le. nlin -nl2 -1) then 

j2=jend IFor bottom mirror! 

!Read new image line: 
line=j -1 

rstatus^rdline (imgin, line, -1 , ibuf f , npixm, arrtyp) 

If (rstatus .ne. SYSNRM) then 

write(6,430) j . . , 

formate "Winbust 042193" :Error reading line ',i4,'.-a 
Call Exit (1) 

Endif 

Do ib=l,nbands ! Overwrite oldest line in barrel-buffer! 
Do i=l,np2 +1 !Left 
i2=i +np2 
i3=np2 +2 -i 

buffn(i2, jend, ib) =ibuf f (i, ib) 
buff n (i3, jend, ib) =ibuf f (i, ib) 

Enddo 

Do i=np2 +2,npix -np2 -1 'Central 
i2=i ■♦'np2 

buffn(i2, jend, ib) =ibuf f (i, ib) 

Enddo 

Do i=npix -np2,npix ! Right 

i2=i +np2 
i3=2*npix +np2 -i 
buf f n {i2, jend, ib) =ibuf f (i, ib) 
buf fn (13, jend, ib) =ibuf f (i, ib) 

Enddo 
Enddo 

Else ! Virtual mirror of bottom of image: 

j2=jend +1 ! Buffer index of "top" line in buffer! 

If (j2 .gt, nl) j2=l 

j2=j2 -1 ! (12-18-92) 

If (j2 .It. 1) j2=xnl ! (12-18-92) 
Do ib-'l,nbands 

Do i»l,npixb 

buf fn (i, jend, ib) -buf fn (i, j2, ib) 

Enddo 



26 



EP 0 624 848 A2 



1^. 



€00 
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800 



Enddo 
Endif 

twrite a status display to standard output! 
if (j -ge. jtype) then 

?±l^ti0x?J ^ine '.i..' i= completed', 

jtype=jtyp€ -rlininc 

Endif 

Enddo iEnd main -j" loop! 

— - call clsimg(iingin) 
Call clsimg (imgout) 
return 
end 
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Job: hist init.vf 

Date: Wed Xpr 21 14:29:59 1993 
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Hist init: (8-12.91) Used in Winbust.042193. Initializes . calcs a his«gra:n 
~ of a window. 

Subroutine Hist init (nlm,np,nl,np2,buf fn,npi^ 
1 npix , nbands , 1 , jend, lb, his , nbinx > 

integer *2 buf fh <npixro, nlro,nbx) 
lnteger*2 his(nbinx) 

integer nlm, np, nl, np2, npixm, npix, nbx, nbands, i, Jend, lb, nbinx 
Integer gmn^gmx, jO 
common /par/gmn,gmx 

...Initialize histogram to zero: 
Do Jt«ginn+l,gmx+l 
liis(k)=0 

Enddo 

Load histogram of window: 

jO^jend 
Do jj-l*nl 

jO-jO +1 

If (jO -gt, nl) jO»l 
Do ii«l,np ^' 

10«=i "np2 -1 +11 

ia«buffn(iO, jO,ib) +1 

his (la) -his (la) +1 

Enddo 

Enddo 

return 
end 
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Job: hist update -vf 

Date: Wed 5pr 21 14:30:52 1993 
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integer«2 buffn (nplwn.nlin,nbx) 

inteier*2 his(nbinx) ^ ,w „hlnx 

integer np.,nl»,np.nl,np2 npi»n,„.x.n.ands,l,ien^ 

integer jO, idrop, iacld,npix 
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jO-jend 

idrop=i -np2 -1 
iadd=i +np2 



il>ffnfidr^oj:^ri^) *1 

fa^^il^nfilid?io:L) .1 
his(ia)-his<ia) +1 

Enddo 



return 
end 
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Job: iselavh 042193. vf 

Date: Wed Apr~21 14:21:43 1993 
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Iselavh: (6-19-91) Mod'd from selavb to eliminate heavy central 21 
weighting on initial-point in interest of speed. 
Should t>e valid for all window sizes. 

Note that this function runs more slowly as the image dynamic range 
[ginax -gmin] increases. 
Assumes histogram is precalculated. 

"ac" is the GL of the "center" pixel in the current window. 
References: 

(1) •Noise Cleaning by Iterattd Local Averaging", by L. Davis 
& A. Rosenfeld, IEEE Transactions Systems, Man, Cybernetics, 
Vol. SMC-8, No. 9, pp. 705-710, Sept. 1978. 

(2) " Digital image smocthing and the sigma filter", 

by J- S. Lee, CVGIP, V. 24, pp255-269 (1983) . Sigma filter, however, 
uses a *fixed* range, in terms of +/- a*sigma, where sigTTia=est. 
std- dev. of additive noise. 

This function computes a selective-averaging filter output for a 
current location of a rectangular window of 1 band of an image array, 
t returns this output as the standard function return. 

Key Parameters: 

his: Int*2 (nbinx) ; Histogram of GLs in current window; 

Assximes that GL value corresponding to a histo bin is 

offset by 1, e.g. his(l) corresponds to GL=0. 
nbinx: Int*4; Max. no. of bins in histogram, 
ac: Int*4; spatially-central GL in current window. 

giiin,gmx:lnt«4; min & max GLs in entire image for present band. 

Passed from caller via Common, 
nav: Int*4; desired no. of pixels to average in the 

current window. 

Output "iselavh" will be a weighted sum of histogram bin occupancies; 
this weighted sum is num/den. 
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Integer Function Iselavh {his, nbinx, ac) 

real* 4 val, vail, val2 

integer*2 his (nbinx) 
integer* 2 cnt 

integer gmin, gmax, nav, ac, nbinx 
Integer ar, al,denr, denl 
INTEGER* 4 NUM,DEN 

common /par /gmin, gmax, nav 

common /par2 /num, den , ar , al , cnt , vail , val2 , val 



50 



! average the nav-nearest grey levels . . . 
den»his (ac +1) 

if (den .ge- nav) then ■ !If have >den local vals =ac,then no-change! 
iselavh=ac 
return 

endif 

num«ac*his (ac +1) 
do )c-l,4999 

ar>ac ! Current GL to "right" of center GL 

al-ac -k !" " "left" of " 

if (ar -le. gmax) then 

num=num +ar*his{ar +1) 
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76 



enddo 



endif 
if (al 



cien=den +his (ar +1) 

. ge. gmin) then 
num=num +al*his Cal +1) 
den =den +his (al +1) 



else if 



22 



endif 
if (den 



.gt. nav) then 
cnt=den -nav 
denr=0 
denl=0 
If (ar 
If (al 



le. 
ge. 



gmax) 
ginin) 



Endif 

write ( 6, ^ 
write (6, ' 
write (6, ' 
write (6, ' 
write ( 6, ' 



denr=his(ar +1) 
denl=his(al +1) 
val2-Float (denr +denl) 
vall=Float (ar*denr) +Float (al*denl) 
val=vall/val2 , , , 

iselavh=Nint ( (Float (num) -cnt*val) /Float (nav) ) 
return 

(den .eq. nav) then 

iselavh=Nint (float (num) /float (den) ) 
return 



Unexpected termination in iselavh!!' 
Local Iselavh values at Exit:' 

Central GL value (ac):',ac 

Gmin , Gmax : ' , gmin , gmax 

Nav: ' , nav 



Call Exit (1) 



end 
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Job: add image. vf 

Date: Wed~Apr 21 14:28:56 1993 
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Add_lmage: (6-26-91) Adds 2 1-3 band images w/ additive 
constants. 

Assumes the internal image representation is Real*4 
Max. # pixels/line=4096. 

Subroutine Add Image (infill, infil2, outf il, f actl, f act2, const, pixformo) 
Character*50 infill, infil2, cutfil 
Real*4 rbuf f 1 (4096, 3) , rbuf f2 (4096, 3) 
ReaX*4 factl, fact 2, const 

integer imginl, imgin2, im9out,nimx, ostatus, gstatus 
integer nbands , nbandsc , pixf orm, pixf ormc, pixformo, band 
integer pstatus, ostatus2, rstatus, wstatus,npx, nln, npx2,nln2 
integer opnimg, getdef , putdef , rdline, wrline 
integer arrtyp, lininc, jtype, imax 

data imax/4096/ 

include ' /local/include/iopackagef ' 
arrtyp=IDREA4 

Open input color image: 

imginl«0 

ostatus=opnimg (imginl , infill, 0, .false. ) 
If (ostatus .ne. sySNRM) then 
write (6, 10) infill 
10 formate In color image ',a50, ' not open-abort') 

Call Exit (1) 

Endif 

gstatus«getdef (imginl, npx, nln, nbands, pixf orm) 
If (gstatus .ne. SYSNRM) then 
write (6, 20) 

20 formate Gstatus error on input — abort') 

Call Exit (1) 

Endif 

If (npx .gt. imax) then 

write (6,*) ' Npx .ge. imax — abort' 

write (6,*)' Npx=',npx 

write (6,*)' lmax= ' , imax 

Call Exit(l) 
Endif 
npix«npx 
nlin=nln 
write (6, 45) 

45 format (/' "Add Image": Npix Nlin Type Nbands Typeo' ) 

write iT, 50) npix, nlin,pixf orm, nbands, pixf oanno 
50 format (13x,2i5, i6, 18, 18) 

— Open 2nd input class map: 
imgin2-0 

ostatus»opnimg(imgin2, infil2, 0, .false.) 
If (ostatus .ne. SYSNRM) then 

write (6, 60) infil2 
60 formate In 2nd image ',a50,' not open-abort') 

Call Exit(l) 

Endif 

gstatus»getdef ( imgin2 , npx2 , nln2, nbandsc, pixf ormc) 
If (gstatus .ne. SYSNRM) then 
write (6, 65) 

65 formate Gstatus error on 2nd input — abort') 

call Exit(l) 
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25 

Call Exit(l) 
Endif 

pstatus-putdef (imgout:,npx,nln,nbands,pxxfonno) 

If (pstatus -ne, SYSNRM) then «„,-r^nr -Abort:' 

write (6,*)' Pstatus error for output-abort 

Call Exit(l) 

ostatus2=opniing ( imgout , out f il , 1 , . f alse . ) 

; "1" for write-only! 

If (ostatus2 .ne, SYSNPM) then ^k«^/ 
ix ^^'^^'^j.i^e^g^*)* Output image not open-abort 

Call Exit (1) 

Endif 

write(6,*)' Factl, Fcat2, Const-' , factl, fact2, const 

Read each line of input image: 
lininc-nlin/4 
jtype=lininc 
Do i»l,nlin 

line=j-l 

•Read multiband line from ?;!p?f® :„„ arrtyp) 

rstatus=rdline (imginl, line, -1, rbuff 1, imax, arrtypi 

\ -"l" for all bands! 
If <r status -ne. SYSNRM) then 

?o^^atr'Er?or reading image line -abort') 
Call Exit(l) 
Endif 

•Read line from I't^and classmap: arrtvp) 
rstatus==rdline (imgin2, line, -1, rbuf f2, imax,arrT;yp; 

I for all bands! • 

If (rstatus .ne. SYSNRM) then 

fo^ltr^lrror reading image line -abort') 
Call Exit(l) 
Endif 

Do band=l , nbands 

^"rb!lffl(i,band)=factl*rbuffl(i,band) 
^ rbuftiUiDanQj^^^^^^^^^^^^2(i,band) +const 

Enddo 

Clip data to output range: 
**if (pixformo .eg. IDBYTE) then 

DO ^'^'"?i2ffi(i,band).Amaxl(rbuffl(i,b^^^^ 

rbuff 1 (i,band) «Aminl (rbuff 1 {l,band) ,255.) 

Enddo , 
Else if (pixformo .eq. IDINTl) then 

DO npix^^^ ^ . ^ ^^^^^ ^) , -128 ) 

rbuff 1 <i,band) -Aminl (rbuffl (irband) ,127.) 

Enddo . 
Else if (pixformo .eq. IDINT2) then 

DO ^"'''^?i$ffi(i,band).Anint (rbuffl (i^ 

rbuffl (i,band) »Amaxl (rbuffl (i,band) , 3^7bo - ) 
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100 



rbuff 1 (i,band)=Aininl (rbuff 1 (i,band) , 32767. ) 26 

Enddo 

Else if (pixfonno .eq. IDREA4) Chen 
Else 

write(6,100) ^ ^ 

formate Pixf or7no=' , 13, ' not impleraented-abort ) 

Call Exit <1) 
Endif 

Enddo !End band loop! 



c Write output line of data: 

line=j-l . 
10 wstatus=Wrline (imgout , line, -1, rbuf f 1, imax, arrtyp) 

If (wstatus -ne. SYSNRM) then 

write (6,*)' Error in wstatus-abort' 
. Call Exit(l) 
Endif 

If (j .ge. jtype) then 

write{6,200) j . , . . ^. ^ 

'5 20O format (1 Ox, ' Line ',i4,' is completed') 

jtype=jtype +lininc 

Endif 

Enddo !£nd line loop! 

z. . .Close images 

Call cisimg (imginl) 
Call cisimg (imgin2) 
Call cisimg (imgout) 
return 
end 
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Bust Class Creates ** two-state ** class map from 3-band residuals image. 
"^Optionally dilates classmap before output. r-i^^^ * ^ 

Calls: Irgb_load_isl, Buster_Class_Create, Dilate_Class, Ioutput_isl. 

Image arrays: , /• 

res : 3-band dirty residuals image (input) 

els :l-band classmap (output) 

residuals image & class map are assumed lnteger*2 values (i.e., not rea 
Note demonstration max image size of 1870 x 1374. 

Does not assume that residuals image is biased by a predetermined 
constant; this info is implicit in the GL-res thresholds entered by the 

Calls: 

Irgb load isl . . ^ • » 

~ (calls image-format -specific lo routines) 
Buster Class_Create 
Dilate'Class 

IoutpuE_isl^^^ image-format-specific io routinesK^ 

Program Bust_Class 

lnteger*2 res (3, 1870, 1374) ^^^^^^ 
lnteger*2 els (1870, 1374) , cls2 (1870, 1374) 

integer*2 npix,nlin 

Integer npixm, nlinm, band, nbands, nbc, nbx, cnt , np, ni 
integer pixf ormi , pixf ormo 

Integer, cnttot 

integer ixO, iyO, nx, ny w,,.-^ 
integer dgl, glmin, glmax, gltmin, gltmax, dilf Ig, byte 

real*4 deltime,Dtime, timearray (2) 

character*50 inf il, outf il, resf il, clsf il 
character*! iname (50) , iname2 (50) 
character* 50 commSO 

equivalence (commSO, iname2 (1) ) 
data dgl/255/ 'Defines "dirt" class! 

data byte/3/ 

data nbx, npixjn,nlinm/3, 1870. 1374/ „^ ^ « 

data pixformo/3/ IClassmap of type "byte 

I — — " 

write (6,*)' [Max. image size: 1870 x 1374)' 
write (6,*)' ' 

format r^Enter in biased residuals image filename [ch503.:',S) 
read (5,' (a50) ' ) resfil 



20 



write (6, 30) ^ , 

30 formate Enter out classmap image name: 

read(5, ' (a50) ' ) clsfil 

write (6,50) ^ w # ci 

50 formate Enter min, max residuals GL thresholds: ,5J 

read (5,*) gltmin,gltmax 

write (6, 60) 
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2 

60 formate Do you wish prefilter dilation? [yes*"!-] : ' , $) 

read (5,*) dilflg 

write (6,*) ' class threshold used is ',dgi 
deltime=Dtime (timeaxray) 

!Read biased color-residuals image: 

Call Irgb__load isl (resf il, res, nJbx,nbands, npixjn, nlinm, 
1 npix,nXin,pixformi, 0) 

write <€, 150) resfil 
150 formate Input residuals file is ',a50) 

I — — ^ 

75 ! Compute class map: 

Call Buster_Class_Creaue (res, cis, nbx, nhands, 
1 npixm, nlinm, npix, niin, gltmin, gltniax, dgl) 

I — 

! Optionally, dilate classmap by one pixel: 

20 If (dilflg .eq. 1) txhen 

write (6,*)' [Dilating class map]' 

Call Dilate_Class (els, cls2, 
1 npixm, nlinm, npix, nlin, dgl) 

Call Ioutput_i5l (clsfii, cls2,npixm, nlinm, 
1 npix, nlin, byte, iname) 

write(6,200) clsfil 
200 formate Out dilated class file is ',a50) 

30 Else ! write undilated classmap! 

Call Ioutput_isl (clsfil , els , npixm, nlinm, 
i npix, nlin, byte, iname) 

write (6, 300) clsfil 
300 formate Out undilated class file is ',a50) 

35 Endif 

deltime=Dtime (timearray) 
write C&, 800) timearray ( 1) , timearray (2) 
800 format (/' Bust_Class:: user time (sec);', 
1 f8,2,'; system time (sec) : ' , f 8 .2/) 
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Job: irgb load isl dugas.vf 
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c lrgb__load_isl: Loads up to a 3-band image into a Int*2 array. 

Subroutine Irgb load_isl (bname,b,nbx, nb, npixin,nlinin, npix, nlin, pixforin, 
1 lunprt)" 
5 Character*50 bnamcinfil 

Integer*2 b (nbx, npixin,nlinm> ,npix, nlin, ibuf f (4096,3) 

integer imgin, gstatus , ostatus 

integer nbands , pixf orm, band, nbx , nb, npixm, nlinm 

integer r status , npx , nln 

integer opnimg, getdef , rdline 

integer arrtyp, lininc, jtype,pixdim 

10 

Data pixdim/4096/ 

include ' /local/include/iopackagef ' 

arrtyp=IDINT2 

c.Open input image: 
IS imgin^O 

inf il=bname 

ostatus«opnimg (imgin, inf il, 0, . false . ) 
If (ostatus .ne. SYSNRM) then 
write (6, 10) infil 
10 formate In image ',a50,' not open-abort') 

stop 

20 Endif 

gstatus*getdef (imgin, npx, nln, nbands, pixform) 
If (gstatus -ne. SYSNRM) then 
write (6,20) 

20 formate Gstatus error on input — abort') 

stop 

25 Endif 

If (nbands *gt. 3 -or. nbands -gt- nbx) then 

write (6,*)' Input no. bands too large — abort' 
stop 

If (npx .gt. npixm .or. npx .gt. pxxdim) then 

write (6,*)' *IRGB_load_isl" :too many pixels—abort' 
30 stop 
Endif 

npix=npx 
nlin=nln 
nb»nbands 
write(6,45) npix, nlin, nbands, pixform 
35 45 format(2x,' "Irgb load isl":: Npix:',i5,'; Nlin:',i5,'; Nb:',i2, 
1 Type:'7i2) " 

If (lunprt .gt- 0) then 

write (lunprt, 45) npix, nlin, nbands, pixform 

Endif 

c Read each line of input image: 

40 linijnc«nlin/4 
jtype«lininc 
Do j«l,nlin 

c Read input lines & load 3-d array: 

line«=j-l 

rstatus^rdline (imgin, line, -1, ibuf f , pixdim, arrtyp) 
! "-1" for all bands! 
45 If (rstatus .ne. SYSNRM) then 

write (7, 80) j 

80 formate Error reading line -abort') 

stop 
Endif 
Do i-«l,npix 

50 
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200 



Do band=l, nbands 

b (band, 1, j) =ibuf f (i,band) 

Enddo 
Enddo 

If (j .ge. jtype) then 
write(6,200) j 

format (lOx, ' Line has been read') 

jtype=jtyp€ +lininc 

Endif 

Enddo 'End line loop! 



c... Close image: 
'0 Call clsimg (imgxn) 

return 
end 
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Job: buster class create_dugas .vf 
Date: Wed Apr 21 17:46:32 1993 
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c Class Createtcreates a fi^-canditote classmap from 
Bie biased color residuals of a "f^|«mooth type aigo^ 
Sir^^eSe^lliLS^irfui^i^^e^^cSnllfe^cfitates of presence 
of defects. 

subroutine Buster ciass.Create (res els, nbx, nbands, npixm. nlinm, 

npix, nlinTgltmin, gltmax, dgl ) 
Integer *2 res (nbx, npixm^nlinm) , els (npixm, nlinm) 
Integer*2 npix,nlin 

integer nbx, nbands, npixm, nlinm, gltmin, gltmax, dgl, glmin, glmax 



75 



10 

11 



25 



fo^M'^sist Class_Create: Gltiain Gltmax Dgl') 



write (6, 11) QfTtmin, gltmax, dgl 
format (19x,2i7,i5/) 

Do j«l,nliii 



Enddo 

return 

end 



Do i-«l,npix . 

glmin-res<l,x,3) 
glmax^'gimin 

Do ib«2,nbands < 41* 

glmin«MinO (glain, res (ib, i, j) ) 
glmax»MaxO (glmax, res (it», i, 3) ) 

??'^tglmin .ge. gltmin .and. glmax .ge. gltmax) then 

cls(i, j)«dgl 

Else ^ ^ 

cls(i, j)=0 

Endif 

Enddo 
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Job: dilate class dugas.vf 
Date: Wed Apr 21 lT:47:31 1993 
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c DiXate Class: Dilates by 1 Pixel a selected class in a 



c 
c 
c 
c 
c 
c 
c 



TO 



! Class: Dilates by i P^^^t t^l^ I -^^^ laulre window- 
ISyte class-map. Dilation kernel is a 3x3 square wi 

ctar«class which is dilated; if current pixel is undilated then 
old classvalue is copied as new classvalue. 

^''|^|!on'Ssing';lrrof'^§xf window still covering image. 

"^roJIt™DnatrcnnT^^ ^^^""^ 
lnteger*2 clsi (npixm, nlinm) , clso (npixm, nlinm) ,npix,nlin 

Integer npixm^nlinm, ctar, im, ip, jm, jp 

I 



10 



T6 



25 



30 



100 



write (6, 10) ctar 

format (/' Dilate_Class: target class iS 

Do j=l,nlin 

jtn=J-l 

jin«MaxO(jm, 1) 

jp-=MinO<jp,nlxn) 
Do i«l,npix 

im=i-l 

iin=MaxO (im, 1) 
ip=i+l 

ip«MinO (ip,npix) 
Do jp 
Do ii»im,ip 



M5/) 



Enddo 



Enddo 



»If any neighbor^ctar, then 

i set center=ctar as well: 

If (clsi(ii,jj) -eq. ctar) then 

clso(i, j)*ctar 

go to 100 

Endif 

Enddo 
Enddo 

!No Change, so copy to output. 
clso(i, j)«clsi(i, j> 
continue 



return 
end 



40 



45 



50 



55 



47 



EP 0 624 848 A2 



10 



75 



20 



999 9 


X XJTX 




y y 


9 99 


ITT* T" 




y y 


9 9 


r 


aaaaa 


y y 


g g 


r 


a a 


y y 


g gg 


r 


a aa 


y yy 


ggg g 


r 


aaaa a 


yyy y 


g 






y 


g g 






y y 


gggg 






yyyy 


i 






t 








t 


li 


oooo 


u u 


ttttt 


i 


o o 


u u 


t 


i 


o o 


M U 


t 


i 


o o 


U U 


t 


i 


o o 


u uu 


t t 


iii 


oooo 


uuu u 


tt 



p ppp 


u u 


pp p 


u u 


p p 


u u 


p p 


u u 


pp p 


u uu 


p ppp 


uuu u 


p 




p 




p 





t 
t 

ttttt 
t 
t 
t 

t t 
tt 



25 



30 







t 




i 






t - 






cccc 


aaaa 


ttttt 


n nnn 


ii 


c c 


a 


t 


zin n 


i 


c 


aaaaa 


t 


n n 


1 


c 


a a 


t 


n n 


i 


c c 


a aa 


t t 


n n 


i 


cccc 


aaaa a 


tt 


n n 


iii 



p ppp 
pp p 
p p 
p p 
pp p 
p ppp 
p 
p 
p 



35 



Job: ioutput isl dugas.vf 
Date: Wed Apr~21 T4:48:35 1993 
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11 

C SUBROUTINE: lOUlirUT isl (10-9-90) 

C PROGRAMMER: ELIZABETH J. DONALDSON/R. Gray 

C DATE: 9/20/90 

C OFILE - INPUT - CHAR*50 - DESTINATION IMAGE 

C FOUT - INPUT - Integer*2 - DESTINATION ARRAY 

C NPIXM - INPUT - INTEGERM r MAX PIXEL DIMENSIONS FOR FIN 

C NLINM - INPUT - INTEGERM - MAX LINE DIMENSION FOR FIN 

C NPIX - OUTPUT - INTEGER* 4 - NUMBER PIXELS IN INFILE 

C NLIM - OUTPUT - INTEGER* 4 - NUMBER LINES IN INFILE 

C PIXFORM - OUTPUT - INTEGER*4 - INFILE DATATYPE 

C INAME50 - OUTPUT^- CHARACTER_-_^UNUSED)^^ 

^^'^ iubroutine"loutput_isl (of ile, f out, npixni,nliiun,npix, nlin^pixf orm, 
+ inameSO) 

include ' /local/include/iopackage . inc' 

character of ile*50, inameSO ( 50) 

15 Integer*2 f out (npixm, nlinm) 

integer*2 npix, nlin 

integer putdef , opnirog, wrband , 

integer*4 imgout , arrtyp, gstatus, pixf orm, npixm, nlinm 

integer*4 ostatus,pstatus, wstatus, npx, nln 

I — — — — 

write (6,5) npix, nlin, pixf orm / ■>» / r.i - ' i4 

20 5 formate loutput isl: image header data: np « ,14, nx - 

+ ' pixf orm = ',i3T 

arrtyp » IDINT2 
imgout=0 

C DEFINE CHARACTERISTICS OF OUTPUT IMAGE 

25 npx=npix 

nln=nlin . 
©status = putdef (imgout , npx, nln, l,pxxf orm) 
if ( gstatus .ne. SYSNRM) go to 1000 

C OPEN OUTPUT IMAGE FOR RANDOM ACCESS . 

ostatus = opnimg (imgout, of ile, 1, .true. ) 
30 if ( ostatus -ne. SYSNRM) then 

write (€,10) of ile 
10 formate Output image' , a50, ' not open. ) 

write (6,*) ' Exiting — ostatus= ostatus 
go to 1000 

endif 

35 ■ if ( pixform.eq.IDBYTE) then 

do 100 j = l,nlin 

do 200 i « l,npix 

if (fout (i, j) .It.O) fout(i,j) = 0 
if (fout (i, j) .gt.255) fout(i,j) « 255 
200 continue 
100 continue 
40 elseif ( pixform.eq.IDINTl) then 

do 110 j -= l,nlin 

do 210 i - l,npix 

if (fout(i, j) .lt.-128) fout(i,j) - -128 
if (fout (i, j) -gt.l27) fout(a,3) " ^.27 
210 continue 
110 continue 
45 elseif ( pixform.eq.IDINT2)then 

do 120 j « l,nlin 

do 220 i • l,npix ^ -,-.,^0 

if ( fout (i,j). It. -32768) /o^t ti, j)--32768 
if ( fout (i, J) .gt. 32767) fout (i, j) ===32767 

50 
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10 



220 continue 
120 continue 

elseif { pixfonn.eq.IDREA4) then 
GO TO 900 

else 

write (6, 30)pixforni 
30 formate lOUTPUT ISL: Pixform ',i3, ' not implemented.') 

endif 

900 vstatus « wrbandCimgout, O,£outrnpixm,nlinin,arrtyp) 
if ( rstatus ,ne. SYSNRM) go to 1000 

1000 continue 

cal 1 cl s img < imgout ) 

return 

end 
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Appendix C — Docket 66,719 
1993 Eastman Kodak Company, Rochester, N-Y. 14650-2201 
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Job: buster erode am dugas.vf 
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25 



- Buster Erode am: (11-11-92) Uses dynamic memory allocation on the Sun. 

: Calls aedirtSxS byte2 nopad; no 2-pixel -dead zone" 

: around image; aXl pixels in image are potentially cleanea. 

: Intended for use on a Sun workstation running SunOS (e.g. a 

: Sparcstation 10) . 

5x5 erosion window is currently hard-wired in. Minimxam clean fraction 
of erosion window is nominally .24, 

Cleans RGB image by grayscale erosion 
using an external 2-state class map. 

Calls: 

Irgb load isl , , ^ 

(Tncludes image- format-specific calls) 
Binput Isl 

" (includes image-format-specific calls) 
Dedi rt 5 x5_byt e2_nopad 
Irgb write isl 

~ (includes image-format-specific calls) . 

Image arrays: 

iml : Original (dirty) 3-band image (input) 

els :l-band classmap (input) 

im2 : Cleaned 3-band image (output) 

Warning: Once erosion iteration begins, iml & im2 contents 
are *partially* -cleaned image; original is not retained. 

Removes by interpolation all pixels whose class-map GL is ,ge. 
a threshold value. Wor)cs on a multiband image (up to 3 bands) . 
Assumes iterative in/of recursive filtering. 

Uses adjustable variable-input window size, allows only filtering 

a rectangular piece of the image. 
Image & class map are assumed Integer*2 values (i.e., not real)^^ ■ 



30 



Program Buster_Erode__Am 

Pointer <piml,iml), (pim2,im2), (pels, els), (pcls2,cls2) 

Integer*2 iml (1, 1, 1) , im2 (1, 1, 1) 
Byte cls(l,l),cls2(l,l) 

lnteger*2 iml (3, 1828, 1332) , im2 (3, 1828, 1332) 
35 Byte els (1828, 1332), cls2 (1828, 1332) 

integer*2 npix,nlin 

Integer npixm,nlinm,band, nbands,nbc,nbx, cnt 
integer pixf ormi,pixformo,pixf ormc 

integer iml_bytes, im2_bytes, cls_bytes, cls2_bytes,tot_bytes 
40 integer itmax,it,cnttot " 

integer inflg 
integer dgl,byte 
Integer*4 iystart, iyend 

real* 4 deltime,Dtime, timearray (2) ,fmin 

45 character* 50 inf il, outf il,clsfil 

data ltmax/40/ !Max, no. of erosion iterations! 

data dgl/255/ ! Defines "dirt" class! 

data diagflg/0/ 

50 
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data byte/3/ 

data nbx,npixin,nlinin/3, 1828, 1332/ 



10 



writl(6;*)' BusterErodelto T'oynaiJiically allocated memory' 

write(6,10) . wca-i.* e\ 

format ( /' Enter in "dirty" file name(ch503 : ,S) 
read (5,' <a50)')infil 



Call Read hdr isl (infil,npix,nlin,nbands,pixf orml) 
write(6,15)npTx,nlin,pixfortni,nbands w:.r,rfc-' 
,0 15 formate npix='.i4,'; nlin-',i4,'; type=',i2,'; bands- ,12) 



npixm=npix 
nlinm^nlin 
nbx=nbands 



15 



20 

20 



30 

25 

40 



30 



write 



iml bytes=2*nbx*npixin*nlinm 

im2~bytes-iml_bytes 

cls~bytes=npixm*nlinn 

tot^5?«f-imp55t^r+im2 bytes +cls_bytes +cls2 bytes 
(6,*) '"Memory allocation requested (bytes): ',tot_bytes 



write <6, 20) , 
formate Enter in class filename:' ,$) 
read (5,' (a50) ' ) clsfil 

formatt'^^^Enter out "cleaned" f ilename [chSO] : ' , S) 
read (5, ' (aSO) ' )outfil 

format e^Enter reqd cleaning fraction "fmin" [Real]:',$) 
read (5,*) fmin 

If (fmin .le- 0. .or. fmm .gt. 0.5) then 

write (6,*)' Fmin out of range — abort 
Call Exit(l) 

Endif 



write'iSj*) ' Class threshold used is ',dgl 

■Read color image to be filtered: 

35 write(6,150)infil . , . , en. 

150 format (/' Input (dirty) file as ',a50) 

deltime=Dtime (timearray) 

plml-malloc (iml bytes) 

If <Pi«^^-f?4^gy.^^^«JJ^ii^, for "iml- failed-abort' 
Call Exit(l) 

Endif 

Call Irgb load Isl (inf il, iml, nbx, nbands,npixn,nlinm, 
1 npix,nTin,pixformi, 0) 

pixformo«pixformi !Out datatype set equal to in type! 

write (6, 160) clsfil 
160 formate Input class file is ',a5G/) 
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pcls=nialloc (els bytes) 
If (pels .eq. oy then 

write(6,*)' Malloc for "els" failed— abort' 

call Exit(l) 

Endif 

Call Binput isl (clsfil,cls,npixm,nlinm,npix,nlin,pixformc, 
l,nEc) 

pi]iiL2"inalloc (im2 bytes) 
If (pim2 .eq, OT then 

write (6,*)' Malloc for *'im2- failed— abort' 

Call Exlt(l) 

Endif 

pcls2=malloc (cls2_bytes) 
If (pcls2 .eq. 0) then 

write (6,*)' Malloc for "0132- failed—abort' 

Call Exit (1) 

Endif 

I 

inflg^l ! Start by assigning "els* is input class map! 

! Clean by erosion iteratively until no more pixels are replaced: 

write (6,*)' [Beginning cleaning]' 
enttot=0 

iystart=^l ! Required initializations 

iyend=nlin 

Do it»l, ittnax 
cnt=0 

If (inflg .eq. 1) then 

Call DedirtSxS byte2_nopad (iml, els, ixn2, 

cls7, npixm, nlinin, npix, nlin, nbx, nbands, 

dgl, ent, fmin) 
Call DedirtSxS byte2_nopad2 (iml, els, iin2, 

1 els"?, npixm, nlininrnpix,nlin, nbx, nbands, 

2 dgl, cnt, fxnin, iystart, iyend) 
lnflg»2 

Else if (inflg .eq. 2) then 

Call Dedirt 5x5_byt e2_nopad ( im2 , els2 , iml , 

els , npixm, nlinxn, npix , nlin , libx, nbands , 

dgl, cnt, fmin) 
Call Dedirt5x5_byte2_nopad2 (im2,cls2,inil, 

1 els, npixm, nlinm, npix, nlin, nbx, nbands, 

2 dgl, cnt, fmin, iystart, iyend) 
inflg^l 

Else 

write < 6,*)' Inflg out of range — atbort' 
Call Exit(l) 

Endif 

cnttot=cnttot +cnt 
write (6,200) it,ent 

formate After it=',i3,' # pels cleaned = ',i6) 

If (cnt .eq- 0) then !No more dirt-class pels to be cleaned! 
go to 500 'Dreaded "go-to** statement! 

Endif 
Enddo 

write($,*)' Warning-Cleaning loop exhausted after it-»' , itmax 
wrlte(6,*)' Aborting due to loop exhaustion' 
Call Exit(l) 
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continue 
write (6, 550) it 

fortnatC Cleaning ended after iteration=' , i.3) 
write (6, 552) cnttot 

formate Total # of points f iltered=' , 16) 

'Release classmap memory: 
Call Free (pels) 
Call Free (pcls2) 

If (Inflg ,eq. 1) then 

Call Free(pim2) Irelease unneeded image memory 

call Irgb__write_isl (outf il, ami, nbXrnbandSr npixm,.nlanm, 

npix,nlin,pixf ormo, 0) 

Else ^ ^ . 

Call Free(piml) Jreleased unneeded image memory 

call Irgb_write_isl (outfil, im2, nbx, nbands, npixm, nlinm, 

npix,nlin,pixf ormo, 0) 

Endif 

deltime=Dtirae (timearray) 
write (6, BOO) timearray (1) , timearray (2) 
format (/' Buster Erode_Am: ; user time (sec):', 

f8.2,'; system time (sec) : ' , f 8 .2/ ) 

write (6, 850) outfil 

formate Cleaned output image file is ',a50) 
end 
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c Binput_isl: (10-7-92) Reads one band of an isl-format image Into 
c a single byte-type array. Input is clipped to range 0-2S5. 



Subroutine Binput isl (inf il, im, npixm. nlinm, npix, nlin, pixf orm, 
band,nban3s) 

Characterise infil 
Character cval 

Byte im <npixrn, nlimn) ,bval 

Integer*2 ibuf f (4096) , npix,nlin 

integer imgin, ©status, gstatus 
integer opniing,getdef , rdline 
integer nbands , pixf orm, band, iband, ival 
integer rstatus, npx, nln, npixm, nlinm 
integer arrtyp, lininc, jtype 



'6 Equivalence (cval,bval) 



c include ' /local/include/iopackagef ' 

! Above comment -out must be removed after f77cvt conversion!!! 

arrtyp=IDINT2 

20 c.Open input image: 

imgin=0 

ostatus»opnimg (imgin, infil, 0, .false. ) 
If (ostatus .ne. SYSNRM) then 

write (6, 10) infil 
10 formate In image ',a50,' not open-abort') 

Call Exit(l) 

25 Endif 

gstatus=getdef (imgin, npx, nln, nbands, pixform) 
If (gstatus .ne. SYSNRM) then 
write(6,20) 

20 format (' Gstatus error on input — abort' ) 

Call Exit(l) 

^0 Endif 

npix=npx 
nlin=nln 
write (6, 45) 

45 format (/' "Binput isl": Npix Nlin Type Nbands Band') 

write (6, STI) npix, nlin, pixform, nbands, band 
50 format (14x, 2i5, 16, iB, iS) 

c.Read each line of input image: 
lininc«nlin/4 
jtypeelininc 
rstatus»0 



^ Do j«l,nlin 

c Read input lines & calc, sum: 

line*j-l 
iband=band -1 

rstatus=rdline (imgin, line, iband, ibuf f , 4096, arrtyp) 
If (rstatus .ne, SYSNKM) then 
45 write (7, 80) j 

80 formate Error reading line -abort') 

Call Exit(l) 
Endif 



Do i»l,npix 
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Enddo 



ival=MaxO (ibuf f (i) , 0) 
ival=MinO (ival, 255) 
cval=Char (ival) 
lm(i, j)=bval 



Enddo 



!End line loop! 



c, -.Close Images 

Call clslmg (imgin) 

return 

end 
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c Dedln:5x5_byte2_nopaa; 9 
c Eliminates "dead zone* 

c of 2 pixels around image edges (i.e. cleaning now takes place at all 

c pixels of the image) . Only need one clean pixel in a 5x5 neighborhood t 

c clean a dirty border pixel (i.e. nroin edge^l) . 

5 ~ 

c Output non-dirt or cleaned classes in classmap have a value of "0". 

c Note that within one call of Dedirt, cleaned GLs are *not* used 

c to clean their neighbors; this will, however, occur in 

c subsequent iterations (further calls to dedirt) . 

c Ad hoc "dirt* remover. Assumes that dirt is detectable as 

70 c a simple grey-level threshold in the input classmap. 

c Filters all pixels *above* a certain 

c threshold grey level ("gth") . 

c Note that a 1-pixel dilation of the classmap is suggested prior 

c to first using this subroutine, to avoid propagation of dirt 

c "ghosts" . 

,5 c 

c Image arrays: (Max dims :npixm, nlinm; Actual dims: npix,nlin) 

c iml: Input (dirty) 3-band image 

c els: Input dirt 1-band classmap 

c cls2: Output dirt 1-band classmap 

c im2: Output (cleaned) 3-band image 

c 

c Input: np,nl rSize of filter window. 

c dgl : Threshold classmap Gh value req'd for cleaning. 

c fmin : Fraction of window needed to be "clean" to be filtered 

c (l=<nmin<np* (nl-1) /2 ). 

c Output :cnt :No. of pixels which were filtered. 

Subroutine DedirtSxS Byte2JNopad (iml, els, im2, cls2,npixm,nlinm,npix,nlin 
25 1 dgl, cnt, fminT 

Integer*2 iml (nbx, npixm, nlinm) , im2 (nbx, npixm, nlinm) 

Integer*2 npix,nlin 

Integer npixm, nlinm, np, nl, np2,nl2, cnt, dent, iO, ix, jO, jx,dgl 
Integer nbx, nb, nxim_neigh, nmin, delx, dely 
30 Integer nmin_edge,nmin_^cent 

Byte els (npixm, nlinm) , cls2 (npixm, nlinm) ,bdgl 

Real*4 sum(3) , wtsum, wt (3, 3) , fmin 
!fmin=.24 implies nrain=6 for 5x5 window. 

35 Character cdgl 

Data np,nl/5,5/ ! Erosion window size! 

Data wt(l,l),wt(2,l),wt(3,l)/l,,l.,0.5/ !Wts specific t: 

Data wt(l,2),wt (2,2),wt (3,2)/l,,0,707,0.4472/ 
Data wt (1, 3) , wt <2, 3) , wt (3, 3) /0 ,5, 0 .4472, 0 . 3536/ 
Data iimin_edge/0/ 



40 



Equivalence (cdgl,bdgl) 



cnt«0 

dcnt-0 

nl2-nl/2 

45 np2«np/2 

nmin_cent«Nint (fmin*np*nl) 
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cdgl^Char (dgl) 10 
write (6, 10) 

10 formate DedirtSxS Byte2_Kopad: Npix Nlin Np Nl Fmin Nmin Dgl' ) 
write (6, 11) npix, nXin, np, nl» fmin, nmin cent, dgl 

11 format <23x,2i5,2i4, f8. 3, i5,i6/) 

! Slide window across image: 
Do j-«l,nlln 
Do l«l,npix 

!i, j are coordinates of center of window! 

!Test local value of classmap: 

If (cls(i,j) .eq. bdgl) then !Filter the point 
dcnt=dcnt +1 
j0=j-nl2 
jx=j+nl2 
jO«MaxO( jO, 1) 
5x«MinO( jx,nlin) 
i0-=i-np2 
ix«i+np2 
iO=MaxO(iO,l) 
ix^MinO (ix, npix) 

!iO,ix,jO,jx are coordinates of extrema of 

! current window 
Do ib=l,nb 

sum (ih) =0 . 

Enddo 

wtsum=0. 

num__neigh=0 

!Calc. over local "dirty" window: 
Do j2=j0, jx 

Do i2=i0,ix 

!If neighbor is *not* dirt, then: 
If (cls(i2,j2) .ne. bdgl) then 

ICalc. a distance-weight for that GL: 

• dis= (Float (i -i2))**2 + (Float (j -j2))**2 

!dis-Sqrt (dis) 

!wt=l,/dis 

! wtsum=wtsxain +wt 

delx»=lnt (Abs(i -i2) +1.) 
dely=Int (Abs(j -j2) +1,) 
wtsum^nrtsiim +vrt (delx, dely) 



Do ib«l,nb 

sum(ib)»siun(ib) +wt (delx, dely) *iml (ib, i. 

Enddo 

! Count no- of non-dirt window pixels: 
num neigh-num neigh -t-l 
Endif - - 

Enddo 
Enddo 

!Calc. local value of "nmin": 

If (j -le. nl2 .or. j .ge. (nlin -nl2 +1)) then 
nmin-nmin edge !Top or bottom border! 

Else if (i .le- np2 .or. i -ge. (npix -np2 +1)) then 
nmin«nmin- edge !Left or right border! 

Else ~ 

nmin^nmin cent !Not in a border region! 
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ndlf 11 

!If have min. no. of non-dirt neighbors, then 
! replace dirty GL w/ wted ave of nondirty GLs: 
If (num neigh .gt. nmin) then ! Replace the center GL: 
5 Do i^l,nb 

ira2 (ib, i, j) =Nint (sum (ib) /wtsum) 

Enddo 

cls2(i,j)=0 'Don't filter point again! 
cnt-cnt +1 

Else !Not enough non-dirt pixels within window! 
10 IKeep the old (dirty) GL for now: 

Do lb=l,nb 

ini2(ib,i, j)«iml(ib,i, j) 

Enddo 

cls2(i, j)«cls{i, j) 
Endif 

75 Else 'Point is non-dirt, no interpolation! 

!Copy input GLs to output: 
Do ib>"l,nb 

im2 (ib, i, j)=iml (ib, i, j) 

Enddo 

cls2(l, j)«cls<i, j) 

. Endif 

20 Enddo !End "i" loop 

Enddo !End - j" loop 

write (6, 500) dcnt^cnt 
500 formate No. of dirt pels:', 16,'; no. of cleaned pels:',i6) 

If (dent .gt. 0 .and. cnt .eq. 0) then 
pc write (6,*)' Warning — stable w/ some dirt remaining!!!' 

write(6,*)'- (aborting]' 
Call ExitCl) 

Endif 

return 

end 
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13 

c li:gb_writ:e_isl : Writes a 3-band image from a Int*2 array. 

Subroutine Irgb_write isl (bname, b, nbx, nb, npixm, nlinm, npix, nlin, 
5 1 pixf orm, lunprD 

Character*50 bname, outf il 

Integer*2 b {nbx, npixm, nlinm) , npix,nlin, ibuff (4096, 3) 
integer inigout,pstatus,ostatus 
integer nbands, pixf orm, nbx, nb, npixm, nlinm 
integer wstatus, npx, nln, ival 
integer arrtyp, lininc, jtype,pixdim, lunprt 
'0 integer putdef, opnimg,wrline 

Data pixdim/4096/ 

include ' /local/include/iopacJcagef 



76 



arrtyp=IDINT2 
nbands=nb 

write(6,45) npix, nlin, nbands, pixf orm , . . . ^ 

45 fonnat<2x,' "Irgb write_isl- ; : Npix:',i5,'; Nlin:',i5,'; Nb:Si2, 

1 Type:'7i2) 
If (lunprt -gt. 0) then 

write (lunprt, 45) npix, nlin, nbands, pixf orm 

20 E^^i^ 

c open output image: 

imgouc=0 
npx=npix 
nin=nlin 
outf il^bname 

25 pstatus=putdef (imgout, npx,nln, nbands, pixf orm) 

If (pstatus -ne. SYSNRM) then 

write (6,*)' Pstatus error for output — abort' 
stop 

Endif ^ . n , 

ostatus=opnimg(imgout,outfil,l, .false.) for write-only! 

If (ostatus .ne. SYSNRM) then 
30 write (6,*)' Output image not open — abort' 

stop 

Endif 

c... Write each line of output image: 
lininc=nlin/5 
jtype==lininc 
35 Do j=l,nlin 

c --. load 2-d array & write lines: 

Do i=l,npix 

Do ib=l, nbands 

ibuff (i,ib)=b(ib,i,j) 

Enddo 
Enddo 

If (pixf orm .eq. IDBYTE) then 
Do ib=l, nbands 
Do i=l,npix 

ival=MaxO (ibuff (i,ib) , 0) 
ibuff (i, ib) *MinO (ival, 255) 

Enddo 
Enddo 

^ Else if (pixf orm .eq. IDINTl) then 

Do ib«l, nbands 
Do i«l,npix 

ival-MaxO (ibuff (i, ib) , -128) 
ibuff (i, ib)«MinO(ival, 127) 

Enddo 
Elnddo 
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20 



100 



14 

Else If (pixrorm .eq. IDINT2) then 
continue 

Else If (pixform .eq. IDREA4) then 
continue 

Else 

write (6, 100) pixform ^ ^ 

formate Pixf onno=' , i3, ' not xmplemented—stop ) 
stop 

Endif 

line-j -1 ... V 

wstatus^wrline (imgout , line, -1, ituf f , pixdxm, arrtyp) 
If (wstatus .ne. SYSNRM) then 

write {6,*)' Error in wstatus — ahort' 
stop 

Endif 

75 c - If (j .ge. jtype) then 

c write(6,200) j ■ ^ 

c 200 format dOx,' Line ',i4,' has been written ) 

c jtyp€=jtype +lininc 



c Endif 

Enddo lEnd line loop! 

c. Close image: 

Call clsimg (imgout) 

return 

end 
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Claims 

1. A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
30 pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

35 d) forming a residual value for each pixel as a function of the difference between the corresponding 

feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
, residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
40 b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

2. A method for the detection and removal of defects in digital images comprising the steps of: 

a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
45 pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

50 d) forming a residuals image which consists of the difference between the feature image pixel values 

and the filtered values of step c); 

e) testing each pixel value of the residuals image of step d) to determine if the residual value is 
within a range of residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
55 b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

3. A method for the detection and removal of defects in digital images comprising the steps of: 
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a) preprocessing a digital image represented by pixel values to form a feature image in which defect 
pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
residuals values expected of defects; 

f) forming a map in which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; and 

g) correcting the digital image using the map created in step f). 

The method according to Claim 1, 2, or 3. wherein the map formed in step f) is conditioned by 
morphological dilation of the pixels marked as defect. 

The method according to Claim 1, 2, or 3, wherein the map formed in step f) is conditioned by 
morphological erosion of the pixels not marked as defect. 

The method according to Claim 1 , 2. or 3, wherein the map of step f) is formed with each pixel marked 
according to the magnitude of its values according to steps b) and e). 

A method for the detection and removal of defects in digital images represented by multiple spatial 
resolutions ordered from lower to higher, comprising the steps of: 

a) preprocessing a lower spatial resolution representation of the digital image by forming a feature 
image in which defect pixels have large local contrast from their neighboring non-defect pixels; 

b) testing the value of each feature image pixel value to determine if the value is within a range of 
pixel values expected of defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects; 

d) forming a residual value for each pixel as a function of the difference between the corresponding 
feature image pixel value and the filtered value of step c); 

e) testing each residual value of step d) to determine if the residual value is within a range of 
residuals values expected of defects; 

f) forming a map In which each pixel that is within the range of expected defects according to steps 
b) and e) is marked as a defect; 

g) correcting the lower resolution representation of the digital image using the map created in step 
f): and 

h) correcting the higher resolution representations of the digital image using the map created in step 
0- 

An apparatus for the detection and removal of defects in digital images comprising: 

means for preprocessing a digital image represented by pixel values to form a feature image in 
which defect pixels have large local contrast from their neighboring non-defect pixels; 

first means for testing the value of each feature image pixel value to determine if the value is within 
a range of pixel values expected of defects; 

filter means for edge-preserving spatial filtering of each of the feature image pixel values; 

means for forming a residual value for each pixel as a function of the difference between the 
corresponding feature image pixel value and the filtered value from said filter means; 

second means for testing each residual value from said forming means to determine if the residual 
value Is within a range of residuals values expected of defects; 

means for forming a map in which each pixel that is within the range of expected defects according 
to said first and said second testing means is marked as a defect; and 

means for correcting the digital image as a function of the formed map. 

An apparatus for the detection and removal of defects in digital images represented by multiple spatial 
resolutions ordered from lower to higher, comprising: 

means for preprocessing a lower spatial resolution representation of the digital image by forming a 
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feature image in which defect pixels have large local contrast from their neighboring non-defect pixels; 

first testing means for testing the value of each feature image pixel value to determine if the value 
is within a range of pixel values expected of defects; 

filler means for edge-preserving spatial filtering of each of the feature image pixel values that are 
5 within the range of expected defects: 

residual means for forming a residual value for each pixel as a function of the difference between 
the corresponding feature image pixel value and the filtered value of said filter means; 

second testing means for testing each residual value from said residual means to determine if the 
residual value is within a range of residuals values expected of defects; 
TO mapping means for forming a map in which each pixel that is within the range of expected defects 

from said first and said second testing means is marked as a defect; 

means for correcting the lower resolution representation of the digital image using the map created 
by said mapping means; and 

means for correcting the higher resolution representations of the digital image using the map 
75 created by said mapping means. 

10. A method for the detection and removal of local defects in digital images comprising the steps of: 

a) generating a feature image from an original digital image; 

b) creating an EPS-residuals image using the feature image; 

c) testing and mapping image pixels based upon the values of the EPS residuals image and the 
pixel values of the feature image; 

d) modifying the mapped image as a function of mapped values of neighboring pixels; 

e) cleaning of defect pixels by an EPS-substitution process or by a greyscale erosion process; and 

f) cleaning of defect pixels in different spatial resolution representations of the original image if 
different spatial resolution representations are present. 

11 A method for the detection and removal of defects in digital images comprising the steps of; 

a) preprocessing a digital image to form a feature image; and 

b) testing each feature image to determine if a value is within a range of pixel values expected of 
30 defects; 

c) edge-preserving spatial filtering of each of the feature image pixel values that are within the range 
of expected defects. 
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